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ABSTRACT 


The objective of this study is to investigate the application of Doppler radar 
systems for global wind measurement. A model of the satellite -based radar wind 
sounder (RAWS) is discussed, and many critical problems in the designing, such as 
the antenna scan pattern, tracking the Doppler shift caused by satellite motion, and 
backscattering of radar signals from different types of clouds, are discussed along 
with their computer simulations. 

In addition, algorithms for measuring mean frequency of radar echoes, such 
as the FFT estimator, the covariance estimator, and the estimators based on 
autoregressive models, are discussed. Monte Carlo computer simulations were used 
to compare the performance of these algorithms. Anti-alias methods are discussed 
for the FFT and the autoregressive methods. 

Several algorithms for reducing radar ambiguity were studied, such as 
random phase coding methods and staggered PRF methods. Computer simulations 
showed that these methods are not applicable to the RAWS because of the broad 
spectral widths of the radar echoes from clouds. A waveform modulation method 
using the concept of spread spectrum and correlation detection was developed to 
solve the radar ambiguity. Radar ambiguity functions were used to analyze the 
effective signal-to-noise ratios for the waveform modulation method. The result 
showed that, with suitable bandwidth product and modulation of the waveform, 
this method can achieve the desired maximum range and maximum frequency of 
the radar system. 
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Chapter 1 


Introduction 


1.0 RADAR REMOTE SENSING 

Since the early 1960s, the field of radar remote sensing has grown into an 
Important technology for scientific research. Radar remote sensing has been 
applied in the areas of 1) astronomical studies, 2) military applications, 3) 
environmental monitoring, and 4) meteorology. The importance of radar remote 
sensing can be attributed to two predominant factors: a) radars do not require the 
sun as a source of illumination, and b) operating in the microwave region, a radar 
signal can penetrate fog, clouds, and to some extent precipitation (rain, or snow). 

During the past, radar remote sensing has found wide applications in 
meteorology, such as storm observation and forecasting. For example, ground- 
based radars are widely used for detecting severe weather and measuring rain-fall 
rate [1-2]. Ground-based Doppler radars, most of them operated in very-high 
frequency (VHF) and ultra-high frequency (UHF), are used for detecting turbulence, 
local wind field etc [3-5]. Recently, S band ground-based Doppler radars were 
developed to detect severe weather, precipitation and velocity fields [6J. As global 
environmental study becomes more and more important, spacebome radar systems 
are anticipated to become more crucial in monitoring global environment, for 
instance, measuring rain-fall rate [7]. 

In this study, we discuss one potential application for radar remote sensing 
— using a spacebome Doppler radar system to monitor the global wind fields. 
Numerous ground-based VHF and UHF Doppler radar systems are being used in 
measuring turbulence and local wind, and recently some microwave radar systems 
were used to measure clouds [8]. However, few spacebome Doppler radar systems 
have been studied for global wind measurement purposes. 


- 1 - 



In this dissertation, we carry out a system study of a spacebome radar wind 
sounder. The major difficulties treated concerning the Implementation of a 
spacebome Doppler radar wind sounder Include power requirements, antenna scan 
pattern, algorithms for measuring mean frequency of Doppler shifts of radar 
signals, removal of frequency ambiguity, and compensating the Doppler shift 
caused by satellite motion. In addition to theoretical analysis, computer 
simulations were used to evaluate the system performance and to compare 
algorithms for measuring mean Doppler frequency of the radar signal and 
algorithms for removal of radar-frequency ambiguity and range ambiguity. 


1.1 BACKGROUND 

1.1.1 NEED FOR GLOBAL WIND DATA 

The importance of investigating a spacebome Doppler radar for wind 
measurement applications arises from the need for global wind data for both 
operational and scientific-research applications. As pointed out in [9], “Knowledge 
of the global wind field is widely recognized as fundamental to advancing our 
understanding and prediction of the total Earth system. Yet, because wind profiles 
are primarily measured by land-based rawinsondes, the oceanic areas (covering 
roughly three quarters of the Earth’s surface) and many regions of the less- 
developed southern hemisphere land areas are poorly observed. The gap between our 
requirements for global wind data and their availability continues to widen. For 
example, as faster computers become available to model the atmosphere with ever 
increasing resolution and sophistication, our ability to do so will be hampered 
because of the lack of data, particularly wind profiles." 

An Improved understanding of the atmospheric wind field is essential for 
purposes such as understanding the physics of the atmosphere, weather forecasting, 
and many others. One of the most important applications for global wind data is in 
numerical weather prediction (NWP), a technique on which modem weather 
forecasts are based. Numerical weather prediction utilizes basic hydrodynamic and 
thermodynamic equations to predict the future states of the environment from the 
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present states [10]. As an Initial value problem. NWP depends critically on the 
accurate specification of the state at time zero. If the models used In weather 
prediction are correct, then the Improvement In weather prediction largely depends 
on the observation and measurement of the state of the weather. This view is 
generally shared by the numerical prediction community [11]. 

Significant progress has been made in NWP in recent years, especially with 
the development of accurate global numerical weather prediction models, improved 
global coverage of the atmosphere provided by satellite observing systems, and with 
the development of high speed computers. However, the current weather forecasts 
are still not close to the theoretical limit of dynamic predictability, generally 
accepted to be about two weeks. Further improvements of weather forecasts are 
considered necessary in the following aspects: the observations that provide the 
initial data for the models, the objective analysis techniques, and the correctness of 
the weather models [12], 

The variables used in NWP are temperature, pressure, and wind. The early 
NWP models were designed to use only pressure and temperature data. Winds were 
derived from the mass observations using the geostrophic relationship. This 
relationship assumes that the latitudinally dependent Coriolis force [13] is 
balanced by the pressure gradient force. This was a natural choice because pressure 
observations were more abundant and more accurate than wind observations. 
Recently, however, it has become increasingly clear that wind data are extremely 
effective for use In numerical weather prediction. Two reasons for this are 
explained by Kalnay, et al. [14]. The first reason is derived from the concept of 
geostrophic adjustment. For most scales of importance to numerical weather 
prediction, the models effectively retain wind data incorporated into the initial 
conditions. Specifically, small-scale pressure-height variations do not result in 
small-scale changes in the wind field. Instead, they are rapidly dispersed as gravity 
waves. In other words, NWP models accept the wind data more readily than mass 
data for scales which can be observed. 

The second reason that winds are an extremely effective source of data comes 
from the well-known fact that integration of noisy data reduces the effect of random 
noise, whereas differentiation enhances the effect of noise. The geostrophic 
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relationship implies that wind is proportional to the horizontal pressure gradient. 
At increasingly smaller scales, the geostrophic relationship is often invalid so that 
wind becomes an increasingly more accurate measure of atmosphere state than do 
the pressure or height measurements. 

In addition, the wind observations are especially Important in the Tropics, 
since the quasi-geostrophic balance, present in the mid-latitudes, breaks down. As 
a result, the wind field cannot be determined from the pressure or height. Moreover, 
a reliable estimate of the divergent component of the wind is necessary to depict the 
convective areas In the Tropics that provide a source of energy for the equatorial 
regions and at times the mid-latitudes. Therefore, wind measurements are more 
important than temperature sounding wherever the winds are not in balance with 
the mass field. This means that they are required to faithfully predict the smaller 
scale systems at all latitudes and all scales in the Tropics. For the larger scales in 
mid-latitudes, temperature data are probably more important, provided they are 
accurate to about ± 1°C. 

Forecast simulations, using wind data in the Tropics and surface wind data 
over the oceans, show significant increases in predictive skill. Global climate 
modeling (GCM) simulation studies show that an rms wind error of 2m/s is 
equivalent to an rms temperature error of about 1°C outside of the Tropics. Thus, a 
wind-measuring system that could achieve such accuracy would be equivalent to the 
best that is possible by any passive temperature measurement system now available 
or under consideration (1 1], 

1.1.2 INSTRUMENTS USED FOR GLOBAL MEASUREMENT OF WIND 

Today's operational wind-velocity observing systems are basically 
Implemented In two forms: those mounted on instrumented towers, and mobile 
Instruments mounted aboard ships, aircraft, or balloons. These instruments are 
very sparse and/ or inherently unable to provide temporal and spatial coverage of 
the global atmosphere at short time intervals. A unique opportunity to measure 
wind with global coverage of the atmosphere is offered by remote-sensing 
instruments mounted on spacecraft. Such instruments could provide observations 
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even In those regions not covered by the conventional ground-based monitoring 
network (examples are the oceans and most of the Southern Hemisphere). 

Efforts have been undertaken by various operational and research centers to 
assess quantitatively, by means of observing-system simulation experiments, the 
potential usefulness of a spacebome global wind sensor and the observational 
requirements that must be met to forecast various atmospheric phenomena [12-14], 
The cost, complexity, and spacecraft constraints indicate that satellite-based wind 
sensors would be most useful In monitoring only the large scale atmospheric 
motions driving the weather systems, while smaller scale motions may be more 
appropriately resolved by ground-based systems. It has been concluded that global, 
twice-daily measurements, with the accuracies and resolutions summarized in 
Table 1.1, would result in more accurate medium-range (up to five days) forecasts In 
the Northern Hemisphere, over which most of the conventional, ground-based wind 
sensors (rawinsonde network) in operation today are concentrated. A major Impact 
is to be expected for forecasting in the sparsely instrumented Southern Hemisphere, 
where the usefulness of forecasts may be advanced by as much as 24 hours [15]. 

Table 1.1 Global and Synoptic Scale Observational Requirements 

Horizontal resolution 100 km (meso-a scale) 

Vertical resolution 1 km (0.5 km in the boundary layer and In 

the vicinity of the jet stream) 

Temporal resolution 6 hour 

Accuracy of the wind component 1-2 m/s in the lower troposphere 

2-5 m/s In the upper troposphere 

Directional accuracy ± 10 degree 


Recognition that only a space-based monitoring system might have the 
capability to provide wind data throughout the troposphere on a global scale has 
prompted researchers to consider new remote wind-measuring techniques. A 
number of active and passive sensing instruments have been proposed, and some 
have already been experimentally tested [16-22], Their main performance 
characteristics are listed In Table 1.2 [15], It is apparent that most of these 
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techniques are unable to satisfy the demands for increasing measurement accuracy 
at better horizontal and vertical resolutions. This Is the case for all passive sensors 
and the scatterometers. On the other hand, excessive size and power consumption 
limit the utilization of Doppler radars for wind measurement from clear air In 
space. 


Table 1.2. Characteristics of various spacebome wind sensors 


Instrument 

Horizon lal 
(km) 

Vertical 

-feTi 

Temporal 



Accuracy 

(m/s) 

Coverage 

Limitations 

Passive 

High-resolution 
Doppler imagers 

125 

4 

24 

< 5 

Middle/upper 

troposphere; 

Stratosphere 

Low reolution and 
accuracy; no tro- 
posphere coverage 

Electro-optical 

modulation 

correlatiors 

150 

5 

24 

< 5 

Stratosphere 

and 

mesosphere 

Low resolution 
and accuracy; no 
troposphere 
coverage 

Cloud-motion 

imagers 

20-50 

none 

0.1 to 1 

2 

10,000 km 2 

No global cover- 
age. No vertial 
profiles 

Active 

Doppler radars 

10 

0.2 

< 0.1 

1 

In precipitat- 
ing systems; 
regional 

Excessive size 
and power 
consumption 

Scatterometers 

25 

none 

12 

± 10% 

Global 

oceanic 

surface 

Low accuracy, no 
vertical profiles 

Doppler lidars 

100 

0.2 to 1 

12 

1 to 3 

Global 

troposphere 

Immature tech- 
nology; no 
coverage in 
cloudy regions 


By comparing the performance of different instruments listed In Table 1.2, it 
can be concluded that Doppler Lldars are viable instruments able to provide, in the 
near future, wind Information as a function of height in the troposphere from clear 
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air regions. The Information could be obtained on a global scale except for the 
cloudy regions, in a grid form suitable for use in operational meteorology, and with 
the required resolutions and accuracies [15]. Since the Doppler radars can also 
provide the required temporal and spatial resolutions, they can fill in the gaps in 
coverage of lidars caused by clouds and rain. Excessive power requirements for 
clear air Doppler radar systems limit the use of radars to wind measurements in 
cloudy regions or precipitation systems. 

1.1.3 REQUIRED OBSERVATION RESOLUTIONS 

Table 1.1 listed the observational requirements for global and synoptic scale 
required in numerical weather prediction modeling. Table 1.3 presents the 
observational resolution generally required to nowcast or forecast various 
mesoscale phenomena, as addressed in a variety of documents (e.g., Federal 
Coordinator for Meteorological Services and Supporting Research, 1982; The 
National Stormscale Operational and R research Meteorology (STORM) 
Program(NCAR, 1984); Shenk et at., 1985; National Environmental Satellite, Data, 
and Information Service, 1985). Typical phenomena associated with the meso-a, 
meso-p, and meso-y scales shown in the table are as follows: meso-a scale 
corresponds to the initiation of a mesoscale convective system(MCS) (NCAR. 1984); 
meso-P scale describes the the internal structure of the MCS (NCAR, 1984), for 
example; and meso-y scale corresponds to small phenomena such as the low-level 
wind shear (microbursts). 

Table 1.3 Required Observational Resolution 


Characteristic Resolution 



meso-a meso-p 

meso-y 

Horizontal 

100 km 

10 km 

0.1 km 

Vertical 

25 mb 

10 mb 

10 mb 

Temporal 

1 hr 

10 min 

1 min 


The last column in Table 1.3, meso-y. is now generally referred to as 
microscale. It follows from Table 1.2 that satellite-based remote sensing platforms 
would be most helpful in meeting the temporal and spatial resolution requirements 
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at the global and synoptic scale and the meso-a scale. For meso-fS and meso-y 
phenomena, ground-based systems such as the Doppler radar wind profilers [11] 
would be more useful. However, in some situations, geosynchronous-satellite cloud 
winds could be useful for meso-JJ events [15], 

1.1.4 THE LAWS FROM EOS 

Since the lidar systems are the most promising remote sensing Instruments 
for global wind measurement in cloud-free areas, a spacebome Doppler laser 
atmospheric wind sounder (LAWS) is under development as part of NASA's Earth 
Observing System (EOS) study [9]. This system is planned for a late-90s launch Into 
low-Earth orbit on an EOS platform. It should make a strong scientific 
contribution to our understanding of the Earth as an integrated system. 

The LAWS will measure the Doppler shift of line-of-sight components from 
aerosol backscattering in the atmosphere with a conically scanned optical 
arrangement. Successive measurements from different directions will provide 
global coverage of wind-vector profiles throughout the troposphere, on a spatial 
scale of 100 km by 100 km at 1 km height intervals, and with an expected accuracy 
better than 1 ms* ^ . Precise management and scheduling of laser pulses should 
allow for more detailed examination of fine-scale meteorological features [23], 

A set of instrument parameters for the lidar system was selected for a Space 
Shuttle orbit of 300 km and a near-polar orbit of 800 km (HufTaker, 1978, 1980). The 
parameters selected for use in the computer simulation are shown in the following 
table. 
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Table 1.4 Base parameters for a Ildar wind sounder 


Space Shuttle Orbit Near-Polar Orbit 


Altitude 

300 km 

830 km 

Target volume 

300 x 300 x 20 km 

300 x 300 x 20 km 

(patch) 



Nadir angle 

62°( 600 km reach) 

52°(1200 km reach) 

Conical scan period 

7s 

19 s 

Wavelength 


9. 1 1 pm ( C 0 2 ) 

Telescope diameter 


1.25 m 

Pulse duration 


6.7 ps 

Optical-detector 


10% 

efficiency 
rms Long-term 


50 prad 

pointing error 
rms Short-term 


2 prad 

pointing error 
Local Oscillator 


50 kHz 


1.1.5 EFFECT ON LIDARS BY CLOUDS 

A large portion of the globe is normally covered by clouds (about 40%), a 
major obstacle for the operation of LAWS. The cloud distribution is usually divided 
Into three layers by meteorologists, (except the vertical clouds: Cumulonimbus or 
Cumulus). The low cloud layer covers the range from the Earth's surface to 2 km. 
These clouds include stratus, stratocumulus, and nimbostratus. The middle layer 
clouds are in the 2 km to 8 km range. These clouds include altostratus and 
altocumulus. The top layer clouds are in the 8 km to 20 km range. These clouds are 
called cirrus, cirrostratus and cirrocumulus. These clouds are generally opaque to 
lidar pulses, except for very thin clouds like cirrus. 
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The research done by HufTaker indicates that 77% of lidar pulses can 
penetrate the upper layer, or down to 8 km [17]. Of those that reach this level, 73% 
will penetrate the middle layer. Thus, only 57% of the transmitted pulses are 
expected to reach down to the 2 km altitude. "Using global cloud statistics, the lidar 
system is not expected to be seriously affected provided that the laser prf is 
sufficiently high. Tropical storms and the warm sector In typical cyclonic storms 
will seriously affect lidar’s performance.” 


1.2 THE DOPPLER RADAR WIND SOUNDER (RAWS) 

The potential for Doppler radar use for global wind measurement 
acquisition was considered poor because of the large antenna aperture and power 
requirements. However, If the Doppler radar is restricted to measuring the wind 
field inside a cloud system, it can provide valuable complementary information for 
the lidar system. In addition, such a system can also provide measurements for 
rainfall rate and ocean surface winds. 

A radar wind sounder (RAWS) is proposed as a multipurpose instrument with 
a scanning pattern similar to LAWS that could measure the winds in the cloudy 
areas [24]. This radar would serve as a complement to the LAWS, which must work 
in clear air or where clouds are thin. Frequencies in X band or Ku band are 
appropriate where precipitation is present, but clouds require higher frequencies 
because cloud drops are very small. More power is needed for clouds than for rain 
alone, but this will be available in future unmanned spacecraft with robust power 
sources. 


1.3 DISSERTATION OUTLINE 

As an Initial step in the study of RAWS, this dissertation is intended to 
investigate several key problems in implementing a Doppler radar wind sounder. 
Several key problems addressed in this study are listed as follows: 
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• Radar backscatter cross section for water clouds and ice clouds. This 
Involves modeling of cloud drop-size distributions and water content, and 
computer simulation of backscattering coefficients from different types 
of clouds. 

• Conceptual design of the system and an analysis of error bound 

* system parameters 

* tracking of satellite speed 

* simulation of signal to noise ratios 

• Algorithms for accurately estimating the moments of the frequency 
spectrum of radar echoes. 

• Algorithms to resolve the frequency ambiguity and range ambiguity 
problems. 

In Chapter 2, we review, as background, the theory of radar scattering from 
particles and clouds, as well as drop-size distributions of water droplets in clouds. 
Computer simulations for a conceptual system are performed to calculate the signal 
to noise ratios for different types of clouds. 

In Chapter 3, we discuss the basic system parameters of a proposed RAWS, as 
well as problems of antenna scan and compensation of the Doppler shift caused by 
satellite motion. In addition, we carry out the analysis of error bounds for 
estimating the wind field for specific situations. 

In Chapter 4, various algorithms for estimating the first moment of the 
Doppler spectrum are discussed. Computer simulations are used to evaluate the 
performance of these algorithms. 

In Chapter 5. radar ambiguity functions are used in analysis of algorithms 
for removal of frequency and range ambiguities. An algorithm using different 
waveform modulations to transmitted pulses is develped to solve the radar 
ambiguity problems of Doppler radars. Again, computer simulations are used in 
this chapter to compare the performance of the algorithms discussed . 
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Chapter 2 


Radar Backscatter and Attenuation From 
Clouds 


2.1 SCATTERING FROM CLOUDS 

In this chapter, we review some of the concepts of radar backscattering from 
particles, as well as some of the models of drop-size distributions of clouds. We use 
computer simulations to estimate signal-to-noise ratios (SNRs) of radar echoes 
scattered from some cloud models. The radar system used In the computer 
simulation is a spacebome Doppler radar wind sounder that is discussed in more 
detail in Chapter 3. Three types of cloud models were selected In the computer 
simulations, according to their mean drop sizes and water contents. The results of 
the simulations are presented as functions of signal-to-noise ratio of radar echo, 
cloud type, and penetration into the clouds. 

2.1.1 DEFINITIONS OF RADAR CROSS SECTIONS 

Assume Si is the power density (W m' 2 ) of an electromagnetic wave incident 
upon a suspended material particle of geometrical cross sectional area A = ?tr 2 . A 
fraction of the incident power is absorbed by the particle, and an additional 
fraction is scattered by the particle in all directions. The ratio of absorbed power P a 
to incident power density Si is known as the absorption cross section 

0,-ff m2 «•» 


Furthermore, the absorption efficiency is defined as the ratio of the absorption 
cross section to the physical cross section of the particle: 
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In addition to absorption, a portion of the incident power Is scattered by the 
particle. The scattering cross-section of the particle, Q s> is defined as the ratio of the 
scattered power P s and the incident power density Si 



and accordingly the scattering efficiency is defined 
Qs 

Ss = 2 (2.4) 

xr 

Both absorption and scattering reduce the Incident power density. The 
extinction cross-section Q e denoting the total power removed by the particle, is 

defined as the sum of the absorption cross-section and the scattering cross section: 
Qe=Qa+ Qs m 2 (2.5) 

The extinction efficiency is defined as 


Se= Sa+ $s (2.6) 

To calculate the power of the radar echoes, knowledge of the backscatter cross- 
section is required. If the back-scattered power density is denoted as Sb. the radar 
backscatter cross-section ob is defined such that the product of ob and the Incident 
power density Sj Is equal to the total power radiated by an equivalent isotopic 
radiator with power density equal to Sb- Therefore, at a distance R from the 
scatterer, the backscatter power density Sb is given as [25]: 



From equation (2.7), the backscatter cross-section of the scatteerr is equal to 

ob = 4 tcR 2 ^ m 2 (2.8) 

b i 

2.1.2 SCATTER FROM A SINGLE SPHERICAL PARTICLE 

Radar cross-sections for objects of almost any shape are difficult to 
calculate. However, the solution for the scattering and absorption of 
electromagnetic waves by a dielectric sphere of arbitrary radius r was derived by 
Mie [26]. The results for the scattering coefficients and absorption coefficients were 
presented in the form of a converging series: 


+oo 

$s(n, X) 4 y/ 2k+ ™ | a k | 2+ | b k I 2) (2 9) 

* k= 0 

and 

^e(n,X)4X (2k + 1)Re(3k+bk) (2 ' 10) 

X k=o 

where x and n are defined as follows 

, 27rr 2;cr i 

x=kbr -^b ~\ 0 y £ rb 

and 



with 
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kb - the wave number In the background medium 
E rb = the real part of the relative dielectric constant of the background 
medium 

Xb = the wave length in the background medium 
Xq = the free space wave length 
np _ the complex index of refraction of the particle 
nb = the complex index of refraction of the background medium 
e C p = the complex dielectric constant of particle 
e c b = the complex dielectric constant of background medium 

The terms ak and bk, known as the Mie coefficients, are functions of n and X- and 
given by Stratton[27) , Battan[28] , and many others. When the background medium 
is air, then e r b = 1 and Xb = Xq- 


2.1.3 RAYLEIGH APPROXIMATION 

When the radius of the particles is considerably smaller than the wavelength 
of the incident wave, specifically when | npx I « 1, the Mie expressions for £ s and 
can be expressed by the Rayleigh approximations (in which only the first and 
second terms of the Mie solutions are considered) (281: 

= (8/3) x 4 I K 1 2 (2.11) 


and 


$ e = 4xIm(-K) + 8/3x 4 I K I 2 + .... 


( 2 . 12 ) 


where 


K = 


n 2 -l 


V 1 


n 2 + 2 e c + 2 


(2.13) 
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n Is the complex Index of refraction of the droplet to the background medium, and £q 
the complex dielectric constant of the droplet relative to the background medium. 
From (2. 12) and (2. 13), the absorption efficiency can be written as 




(2.14) 


The corresponding scattering and absorption cross-sections for a single 
spherical particle are 


Qr ^* 6|K|2 


nr 


(2.15) 


Q a =“X 3 Im(-K) m 2 

The backscatter efficiency in the Rayleigh region is [25] : 


(2.16) 


^ b =4x 4 lKl 2 (2.17) 

Therefore, in the Rayleigh region, the backscattering cross-section of an individual 
spherical particle of radius r is equal to 




64n J 6 


iKh 


m 2 


(2.18) 


The Rayleigh approximations of scattering, absorption, and backscattering 
coefficients are useful In calculating signal- to-noise ratios for radar echoes. The 
only values we need to know are I K| 2 and Im(-K). For water droplets, values of 
I K I and Im (-K) are known for various temperatures and wavelengths: I K I is 

approximately equal to 0.9 for frequencies from 3 GHz to 30 GHz, and temperatures 
from 0°C to 20°C, while Im(-K) increases with frequency [251. Some of the values of 
I K | and Im(-K) for different frequencies and temperature are listed in Table 2.1 

(Table 4.1, Battan[281). 
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The Rayleigh approximation is not always applicable, especially when the 
diameter of the particle is large compared with the wavelength; however, it has been 
shown that the Rayleigh approximation Is valid when I n px I < 0.5. Since, the 
absolute value of refractive Index of a water droplet, ln w l decreases with increasing 
frequency from 1 GHz to 300 GHz [25], the Rayleigh condition In^l <0.5 can be 
satisfied for increasingly larger values of % as frequency increases. The Mie 
extinction and scattering efficiencies, £ e and | s , are shown In Figures 2.1 a, b and c 
as functions of drop radius r [29], The three figures, corresponding to 3, 30, and 300 
GHz, also Include dashed lines which represent the Rayleigh extinction efficiencies. 


Table 2.1 | K | 2 and Im(-K) for Clouds (Table 4.1. Battan [28]) 


Quality 

Temperature 

°C 

10 

Frequency (GHz) 
24.1 

35.5 

1 k i 2 

20 

0.9270 

0.9193 

0.9100 


10 

0.9282 

0.9152 

0.9045 


0 

0.9300 

0.9055 

0.9000 


-8 


0.8902 


Im(-K) 

20 

0.0188 

0.0471 

0.06745 


10 

0.0247 

0.0615 

0.08565 


0 

0.0335 

0.0807 

0.10970 


-8 


0.1036 



The heavy horizontal lines in Figure 2.1 indicate the ranges of drop radii 
characteristic of two types of water clouds and a rain cloud with a rain rate of 25.4 
mm hr‘1. At 3 GHz, the Rayleigh approximation is certainly applicable for the 
water clouds and is approximately valid for the rain cloud, while at 30 GHz, the 
approximation Is valid only for the water clouds. At 300 GHz, the approximation is 
valid only for the fair-weather cloud. 
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a) 



Sphere Rodiuft r (ym) 


b) 

Figure 2.1. Mie efficiency factors for scattering and extinction by a 
water sphere as a function of drop radius (Fraser et al.. 1975 Am. Soc. 
Photogram. [29]). Horizontal arrows indicate drop radii; a) at 3 GHz; 
b) at 30 GHz; c) at 300 GHz. 
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Sphere Rodiui r (pm} 
0 


Ice clouds may contain particles with radii up to about 0.2 mm, but the 
refractive index of ice is smaller than that of water. The combination of these two 
factors leads to the conclusion that, for an ice cloud, the Rayleigh criterion is 
applicable up to about 70 GHz for computing and up to 200 GHz for computing ^5 

1251. 

2.1.4 VOLUME SCATTER 

In a resolution volume in clouds, the scatterers (water droplets or ice 
particles) are assumed to be randomly distributed within the volume such that there 
are no coherent phase relationships between the fields scattered by the individual 
particles. Additionally, the concentration of particles is usually small enough to 
support the assumption that the shadowing of one particle by another may be 
ignored. These two assumptions lead to the conclusion that the total scatter cross- 
section of a given volume is equal to the algebraic sum of the scatter cross-sections 
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of all of the individual particles contained within that volume. Similar statements 
may be made regarding the absorption and backscattering cross sections. 

2.1 .4.1 VOLUME SCATTERING COEFFICIENT 

The volume scattering coefficient k s represents the total scattering cross- 
section per unit volume, and has units of Np m’3 xm^=Np nr*. The volume- 
scattering coefficient k s is given by 

Np m' 1 (2.19) 


r 

max 

v Jp (r) Q s (r) dr 

r min 


where p(r) represents the number of water droplets per-unit volume per increment of 
r, Q s is the scattering cross section for a droplet with radius r. In the Rayleigh 

region, equation (2.19) can be expressed In a summation form 


N 



i=l 


where N is the total concentration of the droplets in a unit volume of the cloud and 
Dj is the diameter of i^ 1 droplet in the unit volume. 

2.1. 4.2 VOLUME ABSORPTION COEFFICIENT 

Similar to the volume scattering coefficient, the volume absorption 
coefficient Is defined as 

r 

max 

K a = Jp(r) Q fl (r) dr Npm’ 1 (2.21) 

r min 

In the Rayleigh region, the absorption coefficient can be expressed as 


- 20 - 


( 2 . 22 ) 


N 

K a = Y Im(-K) Np m _1 

i=l 

Using the following relationship 
N 

m v- 6 *?I D ? S m ' 3 

i=l 

where represents the water content in a cubic meter volume in clouds (the 

fractional volume occupied by the particles multiplied by the density of water( = 10 6 
gm"3), the volume absorption coefficient can be written as: 

k = 6n/X Im(-K) m 10~*> Np m‘ 1 (2.23) 

a V 1 

2. 1.4.3 VOLUME EXTINCTION COEFFICIENT 

For water and ice drops In clouds, the absorption cross-section Q a is much 
larger than the scattering cross-section Q s in the Rayleigh region since Q a is 
proportional to r 3 while Q s is proportional to r 6 . This fact can be observed from 
Figure 2.1. The cloud volume extinction coefficient K e as the sum of k s and K a is 
therefore approximately equal to the volume absorption coefficient x a , and can be 
calculated with equations (2.22) and (2.23). 

2. 1.4. 4 VOLUME BACKSCATTERING COEFFICIENT 

Similar to the definitions of scattering cross-section and absorption cross- 
section, the volume backscattering cross-section is defined as a summation of 
backscattering cross sections of individual drops in a unit volume 
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(2.24) 



= 10-18 



Z 


m 


m 


-1 


where Z is the reflectivity factor (a quantity widely used in the meteorology 

fi _o 

community) In units of mm m" . 

Backscatter from turbulent fluctuations In the refractive index of the 
medium adds to the echo power. However, since the contribution of this type of 
backscatter Is quite small when compared to the droplet backscatter. It can be 
neglected. In addition to the attenuation caused by water droplets in clouds, radio 
waves also suffer attenuation caused by absorption of atmospheric gases. This 
attenuation is mainly caused by the existence of absorption lines of oxygen and 
water vapor in the atmosphere. Oxygen has an isolated absorption line at 118.74 
GHz and a series of close lines between 50 and 70 GHz which act as a continuous 
absorption band. Water vapor has three absorption lines at frequencies of 22.3 GHz, 
183.3 GHz, and 323.8 GHz. As a result, the atmosphere contains a number of 
“windows” where the attenuation of radio waves by atmospheric gases is small. The 
total attenuations caused by absorption of gases as a function of incident angle 
between 20 MHz and 50 GHz are shown In Figure 2.2 [30J. 
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Figure 2.2. Total attenuation of radio waves by atmospheric gases 
versus frequency for various elevation angles (from [30J). 


2.1.5 RADAR EQUATION 

The weather radar equation can be stated as follows [31) 


’r =f C ~2~ P t T ^ ° 2 0 ♦ L M1 11 

Ll024jtMn2 1 1 V J 


expJ 


2 fric, + tc + k ) 
J g P cr 


Ldh 


(2.25) 


where 

c = the speed of light 
Pr = the received signal power 
Pt = the transmitted peak power 

t = the expanded pulse width if the chirp technique Is applied 
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X = the wavelength 
G = the antenna gain 
0 = the horizontal beamwidth 
<|> = the vertical beam width 
Lt = the transmitter loss 
L r = the receiver loss including filters 
r = the distance from the radar to the target 
Kg - the extinction coefficient due to gas 
kp = the extinction coefficient due to rain or snow 
k c = the extinction coefficient due to clouds 

(in the Rayleigh region, it Is equal to k s> the scatter coefficient) 
h = the depth of the radar signal penetration into the cloud 
T| = the volume scatter coefficient 

The total attenuation caused by the atmospheric gases is dependent upon the 
pointing angle of the antenna and the operating frequency. At a pointing angle of 30 
- 35 degrees and frequency of 35 GHz, the total loss is less than 0.5 dB, as shown in 
Figure 2.2. Kp and k c are dependent upon cloud type and can be calculated using 

equation (2.22) or (2.23). 


2.2 MODELS OF DROP-SIZE DISTRIBUTION IN CLOUDS 
2.2.1 CONTENTS OF CLOUDS 

Cloud droplets are usually formed by water vapor condensing on particulate 
which serve as condensation nuclei. Supersaturation (humidity of >100%) is 
required to condense water vapor In pure air. For a visible cloud, the droplet's 
diameter is > 5 micrometers. The concentration of water droplets in clouds is on the 
order of 100 per cubic centimeter and typical radii are about 10 pm. This structure is 
extremely stable as a rule, and the droplets show little tendency to come together or 
to change their sizes except by general growth of the whole population [32]. In 
addition to mean radius of droplets, total concentration and water content 
per-cubic meter are other important parameters used to classify cloud. s In 
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non-precipitating cumulus, a typical value for water content is 0.59 gm-3, with peak 
values of about lg m'3. in stratus clouds, the values tend to be smaller. In 
cumulonimbus clouds, the water contents can exceed 5 gm'3. 

In addition to water droplets, clouds often contain ice particles when the 
temperature Is low. However, observations show that the 0°C level in the 
atmosphere does not induce a sharp discontinuity in the microstructure of clouds 
and precipitation. Clouds have a high probability of containing no ice if the 
temperature of the cloud is warmer than -10°C. However, with decreasing 
temperature the likelihood of ice increases: below -20°C, more than 90 percent of the 
clouds contain ice particles [33], 

2.2.2 DROPLET SIZE DISTRIBUTION OF WATER CLOUD 

Many cloud physicists have published results of measurements of drop-size 
distributions or liquid-water contents, or both, in various types of cloud and fog. 
Some of these results are shown in Figures 2.3, 2.4, and 2.5. These examples show 
that most drop-size distributions measured in many different types of cloud under a 
variety of meteorological conditions exhibit a characteristic shape. Generally, the 
concentration rises sharply from a low value to a maximum, and then decreases 
gently toward larger sizes; thus, the distribution becomes positively skewed with a 
long tail toward the larger sizes. Such a characteristic can be approximated 
reasonably well by some analytical formulas [34-35]: 
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1.0 



Droplet Diameter (pm) 


Figure 2.3. Measured drop-diameter histograms for different kinds of clouds. 
Note the change in ordinate scale from part to part. Note significant numbers 
of large drops in all but (d). (a) Orographic cloud. Hawaii, p = 0.4 gm-3. (b) Dark 
stratus over Hilo, HI, p = 0.34 gm-3. (c) Tradewind cumulus over Pacific near 
Hawaii, 615 m above base, p = 0.5 gm-3. (d) Continental cumulus over Blue 
Mountains near Sydney, Australia, 615 m above base, p = 0.35 gm-3 (From [35]) 



Figure 2.4. Average cloud-drop spectra reported by aufm Kampe and 
Welckmann for different cloud types. Note the large number of large drops 
present in the cumulus congestus and cumulonimbus clouds (from (351). 
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Figure 2.5. Distribution of drop radii for summer convective clouds 
over the U.S. based on work of Battan and Reitan. Left: 19 fair- 
weather cumulus clouds with average of 293 drops cm-3. Right: 
cumulus congestus for two cases: arrested growth with 247 drops cm-3 
and those growing to point where 1950s weather radar showed echoes 
with 188 drops cm-3 (from [35]). 


LOG normal distribution 


n(a) = ■ 


f/2 


■sxp) - 


ln*(a/a m ) 


jta 


2a 2 


(2.26) 


where a is drop radius, am is the median drop radius and a is the root mean square 
(RMS) deviation of the logarithm of the drop radius. 

MODIFIED GAMMA DISTRIBUTION 

n(a) = A aP e' Ba7 (2.27) 


where A, B. p. y are positive parameters. The maximum of this distribution occurs 
at am. an observable quantity which relates the parameters of the distribution to 
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A constraint is imposed on equation (2.27) in terms of the total concentration N. 


AB -<P +1 W 

r«p + d/y) 


(2.28) 


KHRIGIAN AND MAZIN DISTRIBUTION 

The Khrigian-Mazin drop-size distribution is a special case of the modified 
Gamma distribution [34], which can be expressed as 

n(a)=Aa 2 e(- Ba ) (2.29) 


The parameters A and B can be related to the first and second moments of the 
distribution: 


N = 


r max 
Jn(a) da = 
r min 



(2.30) 


and 


rmax 

<a> = 1/N Jan(a)da = 3/B (2.31) 

main 

where <a> is the mean radius of drop size. Another related quantity of interest is 
the water content, Wl. For the Khrigian-Mazin distribution, 

A = 1.45xl0 18 -^ k ^ (2.32) 

p <a>° 

wl 

p <a> 8 


N= 1.07 x10 s 


Number of Drops 


(2.33) 


W L = 0.934579 x 10' 5 x N p <a> 3 g m' 3 


(2.34) 


where <a> Is In pm, p Is in g cm' 3 , Wl Is in gnr 3 and N Is In number of drops per cm 3 . 

The Khriglan-Mazln distribution Is very convenient to be used In computer 
simulation of cloud drop size distributions and calculation of water contents of 
clouds while the modified Gamma distribution may give overflow problems in 
simulating drop-size distributions or some types of cloud. The log-normal 
distribution does not have the simple relationships between the total concentration 
of droplets, the water contents, and the mean radius like the other distributions. 
However, all of the analytical expressions given represent only average 
distributions. Individual drop-size spectra may be significantly different. Figure 
2.6 shows computer simulated Khrigian-Mazin Drop size distributions for three 
different types of clouds. These three models of clouds are classified as thin, 
medium, and heavy according to their water contents per-unit volume. 



Figure 2.6. Computer simulated drop-slze-distrlbutlons of three different 
types of clouds 0.3, 0.5, and 1.0 gm' 3 (thin, medium and heavy). 
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In general, when the drop size Increases, the concentration of droplets 
decreases. The following parameters may be helpful in cloud modeling [32]: 


Table 2.2. Some Typical Parameters for Clouds 


Droplet Type 

Mean Radius 
(pm) 

Concentration of 
Droplets (per liter) 

Falling Velocity 
cm s' 1 

Typical Cloud Drop 

10 

10 6 

1 

Large Cloud Drop 

50 

10 3 

27 

Typical Rain Drop 

1000 

1 

650 


2.2.3 ICE CLOUD MODELING 

The drop size distributions of ice crystals In clouds are not as well known as 
those for water droplets. With the average size of the ice crystals under many 
conditions being considerably large and their shapes irregular and usually far from 
spherical, it is difficult to model an ice cloud. Ice crystals in ice clouds can attain 
sizes an order of magnitude larger than water droplets found In water clouds. Hence, 
the reflectivity factor Zi of an ice cloud may be several orders of magnitude larger 
than that of a water cloud with the same liquid water content %. 
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Figure 2.7. Average ice crystal spectra in a) Ci Spi. b) AS, c) and d) Ci 
unc, e) Cs, f) Ac and g) Cb cap. The size class is 200 p (From [36]). 
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The water content of a cloud Is typically less than 1 gm-3 and rarely exceeds 
4 gm-3. The factor | K w 1 2 for water varies between 0.89 and 0.93 over the 0-20 °C 
temperature range and 1-10 cm wavelength range. For ice, |Ki| 2 is about 0.2, 
which is 4.5 times smaller than | K w | 2 . but because of the much larger Zt ( in 
comparison to Zy,), ice clouds are much more readily detectable by radar than water 
clouds. 


Figure 2.7 shows some averaged ice crystal spectra in different cloud types 
measured by Heymsfield and Knowllenberg [361 . For cirrus clouds, Heymsfield and 
Knollenberg measured the following average characteristics: 

• ice crystal concentration 1.0 xlO 4 to 2.5 x 10 4 m® 

• mean crystal length 6.0 xlO' 2 to 1.0 xlO’ 1 cm 

_Q 

• ice-water content 0. 15 - 0.25 gm 

• radar reflectivity factor 5.0 - 20.0 mm® m ® 

• precipitation rate 0.5 - 0.7 mm hr* 1 

2.3 COMPUTER SIMULATIONS AND CONCLUSIONS 

Most of the early experiments used 10 cm to 3 cm radars that could only 
detect drops larger than a few hundred microns in diameter. More recent high- 
power 3-cm radars, and most 1-cm radars, permit the detection of drops with 
diameters larger than a few tens of microns. A majority of this work has shown 
that the appearance of the radar echo is characteristically related to the cloud 
dimensions and temperatures [37-38], There is a correlation between the number of 
echoes and the cloud dimension and top temperature. Larger clouds produced more 
radar echoes. Similarly, louds with colder cloud top temperatures also produced 
more radar echoes. The explanation for this phenomena is that the large clouds 
tend to have a broader spectra of drop sizes and greater Z factors with the same 
mean drop size and water content. Cold cloud temperatures indicate that the clouds 
may contain ice particles and therefore have large reflectivity factors. 
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2.3.1 COMPUTER SIMULATIONS OF SNR OF RADAR ECHOES 


Three cloud types, rated as thin, medium and heavy according to their water 
contents and mean drop sizes, have been selected In the computer simulations. 
Their water contents and mean drop-radii are listed In Table 2.3. The basic purpose 
of these computer simulations was to examine the signal to noise ratio of radar 
echoes from different cloud types as a function of cloud penetration, frequency and 
reflectivity factor. 


Table 2.3. Parameters In Cloud Modeling 


Type of Cloud 

Water Content 
gm' 3 

Mean Drop Radius 
pm 

Thin Cloud 

0.3 

4.5 

Medium Cloud 

0.5 

7.5 

Heavy Cloud 

1.0 

11.6 


A menu-driven simulator was developed, which allows the user to set up 


parameters, such as frequency, antenna gain, altitude, signal-to-noise ratio, cloud 
models, etc., from popup menus. The parameters for the radar system used in the 
simulations are listed in Table 2.4. A detailed discussion of these parameters is 
presented in Chapter 3. 

Table 2.4. Parameters Used in Computer Simulation 

Parameter 

35 GHz 

10 GHz 

Antenna Gain 

68 dB 

57 dB 

Beamwidth 

0.00122 rad 

0.00427 rad 

Peak Power 

3000 W 

3000 W 

Pulse Length (chirped) 

1 ps 

1 ps 

Chirp Gain 

20 

20 

Pointing Angle 

30 deg 

30 deg 

Transmitter Loss 

1.5 dB 

1.5 dB 

Receiver Loss 

1.5 dB 

1.5 dB 

Noise Figure 

4 dB 

4 dB 

Im(-K) 

0.08565 

0.0247 

K 2 

0.9 

0.9 
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using the parameters given in Table 2.4. computer simulations of slgnal-to- 
notse ratio (SNR) versus cloud depth were performed for different types of clouds 
requencles. and altitudes of orbits. The mean radius of droplets, drop 
concentration. and water contents of a cloud are no. Independent. An example of 
e corre atlon between these values lshown In Figure 2.8 1351. From this figure an 

experimental formula was derived that relates the mean drop sfee and the water 
content In a cloud. 


r = 11.6 + 13.5* log(W) 


(2.39) 


where r ls mean radius, and W Is water content. This equation Is only used for the 
purpose of computer simulations of cloud drop-slxe distribution hr this chapter 
The actual clouds may have a much different relationship than (2.39). 





The results presented tn Figures 2.9 to 2. 1 1 show that, for water clouds the 
return signal from cloud tops a. 35 GHz ls much larger than a. 10-GHz. This Is due 

°' aC , < ' s' thC 35 GHZ SySten> has a h ‘ ghCT anl ' nna « aln “<1 larger backscatter 
coefficients than the 10 GHz system does. 
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Figure 2.9 shows the slgnal-to-nolse ratio for the “thin" cloud type (like fair- 
weather cumulus). In this condition, the radar echo is weak even at 35 GHz: 5 dB 
SNR at a 300-km orbit, and below zero dB at a 800 km orbit. For the 10-GHz system, 
the signal is too weak to be useful in measurement of wind vectors from clouds. 
Figure 2.10 shows the slgnal-to-nolse ratios from a medium cloud, with water 
content 0.5 gm"^.; the stgnal-to-nolse ratios increase substantially (about 10 dB) as 
the water content increases from 0.3 gm'3 to 0.5 gm'0-3. in this case, at both the 300- 
km and the 830-km orbits, the 35-GHz system is able to provide high enough SNRs 
for measuring Doppler frequencies. As for 0.3 gm‘3, the 10-GHz system cannot 
provide high enough SNR Figure 2.11 shows the signal return from a heavy cloud 
with water content 1.0 gm‘3. The SNR’s at both frequencies are further increased: 
even at 10 GHz the SNR is above 0 dB for the 300-km orbit. 

From these results, it can be concluded that a frequency of 35 GHz or higher 
is required for measuring moments of Doppler spectra of radar echoes from a 
majority of water clouds, as well as ice clouds. Although the 10-GHz system shows 
greater penetrations than the 35-GHz system, especially in the heavy cloud (1.0 
gm‘3), it is limited for wind measurement In clouds from space because of low SNR 
and wide antenna beamwidth (as we will see in Chapter 3). A frequency around 10 
GHz can be used for measuring wind in rain. The rain may cause too much 
attenuation for 35 GHz or higher. In addition, the 10-GHz radar can also determine 
ocean surface wind fields. 
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Figure 2.9. Signal to noise ratio of radar echo as a function of cloud 

penetration (water content of cloud = 0.3 gm' 3 ); a) orbit height = 300 
km; b) orbit height = 830 km. 
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Figure 2. 10. Signal to noise ratio of radar echo as a function of cloud 
penetration (water content of cloud = 0.5 gm"3); a) orbit height = 300 
km; b) orbit height = 830 km. 
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Figure 2.1 1. Signal to noise ratio of radar echo as a function of cloud 
penetration (water content of cloud =1.0 gnr 3 ); a) orbit height = 300 
km; b) orbit height = 830 km. 
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Figure 2.12. Signal to noise ratio of radar echo as a function of 
reflectivity factor Z (mm® m"®). The return signal is from 1 km depth 
from the cloud top. 


Finally, Figure 2.12 shows the simulated signal-to-noise ratio (SNR) versus 
increments of the Z factor. Because the ice cloud is difficult to model, we can derive 

some qualitative idea from Figure 2.12 about the signal to noise ratio from ice 
clouds. For example, although the | Kj | is about 0.2, 4.5 times smaller than 

| K w | , Z for cirrus (ice) clouds at 35 GHz is 5 to 20 mm 6 m'®, and the signal-to- 
noise ratio is above 20 dB based on Figure 2.12. 

2.3.2 CONCLUSION AND FUTURE WORK 

The results presented in this chapter demonstrate that, from a SNR point of 
view, the 35 GHz Doppler radar can provide high enough SNR’s of radar echoes from 
clouds for measuring mean Doppler frequency. Although the computer simulations 
were based on water cloud models, the results may also be applicable to ice clouds as 
discussed in the previous section. However, many other factors were not considered 
in the computer simulations and therefore not discussed in this chapter. For 
example, we did not consider the cases that the cloud dimensions are smaller than 
the antenna beamwidth, or the thickness of clouds are narrower than the vertical 
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resolution of the radar. In these cases, the signal-to-nolse ratios of radar echoes 
would be smaller than those simulated in this chapter. Moreover, some of the cloud 
models used in the computer simulations may be Impractical. These questions need 
to be addressed in future studies. 
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Chapter 3 


Conceptual Design and Analysis of 
Performance 


3.1 INTRODUCTION 

A large volume of published literature concerns operations of Doppler 
radars for measuring wind speed, weather forecasting, detecting turbulence, etc. 
However, there was not conclusive evidence that Doppler radars were able to 
measure the wind fields in clouds from space. One major obstacle was that the 
radars lacked the power to overcome the weak backscatter from clouds. Although 
recently there has been research in the area of using VHF and UHF Doppler radars to 
measure wind speeds [39-42], these VHF and UHF systems could not be used in 
spacebome applications because of the excessive size of the antenna and power 
requirements. To meet the requirements for vertical resolution and peak power, a 
spacebome Doppler radar for wind measurement should use a high frequency, as 
discussed in Chapter 2. 

The Radar Wind Sensor (RAWS) Is a proposed space-bome Doppler radar to 
be used as a multi-purpose Instrument for measuring global wind fields in cloud- 
covered areas as well as rainfall and surface winds on the oceans. The basic 
configuration, as shown in Figure 3.1, has two antenna beams at two fixed elevation 
angles 4>1 and <J>2. and the antenna scan patterns are conical. The transmitted pulses 
are frequency modulated with compressed pulse widths of about 1 (is. Each antenna 
beam transmits signals at two different frequencies: one frequency is at about 10 
GHz (X band), and the other is approximately at 35 GHz (Ka band). The frequencies 
used on the two antenna beams may need to be slightly different to avoid 
interference between received signals. To concentrate on the topic of wind 
measurement from non-precipitating clouds, we will primarily discuss the system 
configuration and performance for 35 GHz. The 10-GHz system will be covered by 
others since It is primarily useful for measurements of wind in precipitation 
systems, ocean-surface wind, or rainfall rate. For the lower elevations in rain, a f 
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requency like 10 GHz is required because of the high attenuation of the 35-GHz 
signal traveling through rain. 



Figure 3.1. Basic concepts of two beam conical scan: the same area 
can be viewed forward and rearward by the antenna with different 
looking angles for wind vector component retrieval. 
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With circular scanning patterns, areas of hydrometeor activity can be 
illuminated multiple times by the antenna beams when the satellite traverses over 
these areas. Since wind vectors are three dimensional, one needs at least three 
measurements of the Doppler shills from different angles to retrieve a wind vector. 
In clear air or light clouds, it is often assumed that the vertical components of wind 
vectors are uncorrelated with the horizontal components, and that the amplitudes 
of vertical components of wind vectors are very small in comparison to the 
horizontal components. Under such assumptions, the vertical components of wind 
fields can be estimated separately from the horizontal components, or even Ignored 
In calculation of the horizontal components of wind vectors. Therefore, it is 
possible, by measuring the Doppler frequencies of water droplets in these regions 
from two different angles, that the horizontal wind vectors in the resolution volume 
can be retrieved. However, in heavy cloud regions and precipitation systems, the 
vertical wind vectors may not be ignored. In such cases, measurements with three 
linearly Independent pointing angles are necessary to retrieve the wind vectors. 

An assumption concerning the wind field was mentioned in Chapter one: the 
RAWS is used to measure the global and meso-a scale phenomena of the 
atmosphere. Thus, the measurements made within a 100 km by 100 km area can be 
averaged to achieve a more reliable measurement of the wind vector in this area. 

In Chapter 2, we addressed the radar backscatter from different types of 
clouds. In this chapter, we will discuss the system parameters of RAWS and conduct 
an error analysis of the system performance. In addition, we will also discuss 
strategies of antenna scan patterns, and a method for compensating Doppler shifts 
caused by satellite motion. 

3.2 BASIC CONFIGURATION OF THE SYSTEM 

The basic parameters of a conceptual system for the radar wind sounder are 
listed in Table 3.1. Tradeoffs encountered in selecting the parameters such as 
antenna elevation pointing angle, scan period, and scan modes are discussed In 
detail below. Two orbit heights were selected for the conceptual design of the radar 
wind sounder. One is at 300 km for space-shuttle orbit, and the other is at 830 km 
for near-polar orbit. These orbit selections were referred to in [9]. In addition, some 


- 42 - 



major obstacles, such as clutter rejection and limit on vertical resolution, 
identified in the following discussion. 


are 


Table 3. 1 Basic Parameters for the Radar Wind Sounder 


PARAMETER 

10 GHz 

35 GHz 

Altitude of Satellite 


300 km or 830 km 

Target Volume Used In Output 100 km x 100km x 20 km 

(from many Individual measurements) 

Nadir Angle 30° and 35° 

PRF 


3500 Hz 

Pulse Width (Compressed) 


1 yis 

Time-Bandwidth Product 


20 

Antenna Size (Parabolic) 


8m 

Antenna Gain 

57 dB 

68 dB 

Horizontal Beamwldth 

0.00427 rad 

0.00122 rad 

Vertical Beamwidth 

0.00427 rad 

0.00122 rad 

Scan Period 

10 s 

10 s 

Vertical Resolution 

2 km 

1 km 

Footprint at 30° (300 km) 

1.5x 1.7km 

0.5 x 0.4 km 

Footprint at 30° (830 km) 

2 x 5. 1 km 

1.5x 1.2 km 

Doppler Bandwidth due to 
Satellite Motion 


1900 Hz 

Peak Power 


3000 W 

Average Power 


210W 

Receiver Noise Figure 


4dB 

Transmitter Loss 


1.5 dB 

Receiver Loss 


1.5 dB 

Spacecraft Velocity 


7.5 kms'l 
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3.2.1 SYSTEM PARAMETERS 


OPERATING FREQUENCY 

The operating frequencies of RAWS were chosen to be at 10 GHz and 35.5 GHz. 
However, from the computer simulation study In Chapter 2, we concluded that 
operating frequencies around 10 GHz may be more appropriate for measuring winds 
in precipitation, surface winds over the ocean, and rainfall rate. To measure winds 
from clouds, frequencies around 35 GHz or higher are needed with current 
technology. By observing the transmission windows in the microwave region 
presented in Figure 3.2, for frequencies above 20 GHz, windows at frequencies 35.5 
GHz, 90 GHz, or 135 GHz are all applicable for measuring wind fields in clouds. 
Higher frequencies can increase the SNR of radar echoes from thin clouds and 
reduce the antenna size, but they also endure more extinction, and may not be able 
to penetrate deep enough through heavy clouds. In this dissertation we will only 
discuss the case of 35 GHz. The feasibility of using higher frequencies to measure 
wind fields from space may need further investigation. 

ANTENNA GAIN, SIZE, AND BEAM WIDTH 

To provide the high gain level and narrow vertical beamwidth required for 
good vertical resolution, a large antenna is desired. With the technology currently 
available, the size of an antenna can be made as large as 1000 wavelengths. We 
chose to use an 8 m diameter parabolic antenna in this conceptual design. It Is 
possible to use this antenna for both 35 GHz and 10 GHz frequencies. For a 
uniformly illuminated parabolic antenna, the gains of the antenna for both of these 
frequencies can be calculated with the following formula: 

G 0= D O e ap e t 

= (t) e ap e t (31) 

where Go is the maximum value of gain. Do is the directivity, d is the diameter of the 
dish of the reflector, eap is the aperture efficiency of the parabolic antenna, and et is 
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the antenna efficiency. Typical values for these parameters are e a p = 60% - 80% 
andet=95% [43]. If we choose e a p=70%, et =95%, and substitute these values In to 
(3.1), It follows the antenna gains of G=68 dB at f=35.5 GHz and G = 57 dB at f=10 GHz. 


WoveUngth (cm) 



Ulaby. R K Moore, and A. K. Fung. 1981 (25D 


The beam widths of the antennam can be calculated using the experimental 
formula from (43): 


g Q = 


30000 
®ld ®2d 


(3.2) 


where 0id and 02d are horizontal and vertical beam widths. For a parabolic 
antenna, from (3.2), 0 1 d = 02d = O.O7°=1.22 m rad at 35.5 GHz and 
01d=02d=O-244°=4.27 m rad at 10 GHz. 
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ANTENNA POINTING ANGLE 


There are tradeoffs in selecting the antenna pointing angles, the vertical 
resolution and the swath width. The vertical resolution imposes a limit on the 
use of large antenna pointing angles; however, a large pointing angle is desired 
for achieving large swath width and reducing measurement errors for 
horizontal components of wind vectors. Figure 3.3 shows the limitation on the 



Figure 3.3. Tradeoff between the antenna pointing angle, antenna 
beamwidths and the vertical resolution, where <|> is antenna pointing 
angle, r e is the edge-to-edge vertical resolution, and H is the height of 
the antenna. 

antenna pointing angle imposed by the vertical resolution. The vertical resolution 
not only depends on the pulse length, but also on the antenna pointing angle and 
beamwidth. Assume that R is the slant range, <|>i the elevation angle of the antenna. 
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3 the beam width, x the pulse length (compressed, if chirp is applied), then r e (the 
edge-to-edge vertical resolution) is equal to: 

r e = 3 R sin y cos i = 1,2 (3.3) 


Equation (3.3) shows that the vertical resolution is a sum of two terms. The first 
term is contributed by the antenna beamwidth, and the second term is contributed 
by the pulsewidth. For the assumed parameters : 


R = 830 km /cosft 
r e = 1 km 

3 = 1.22 mrad (35 GHz) 
4>1 = 35° 


the pulsewidth is limited by the following inequality: 


x < 2 


r e - 3 R sin 4>. 
c cos <|>. 


(3.4) 


From (3.4), the pulse length x is limited to less than 1.65 ps. Using the same 
parameters as above except an orbit height cf 300 km (space-shuttle oihit), the pulse 
length is limited to less than 4.08 ps. However, for f=10 GHz and 3=4.27 mrad, the 
vertical resolution cannot be less than 1 km in 830-km orbit and the pulse length is 
limited to less than 0.56 ps (to satisfy 1 km vertical resolution in 300 km orbit). 
The primary limit to resolution is the beam width. Hence, a taller antenna (>8 m) 
would allow improved vertical resolution. 

Another definition of the vertical resolution is called the effective vertical 
resolution. It is defined as the -6 dB 3-dimensional contour of the radar 
illumination [44], With this definition, Kozu (1989) pointed out that the effective 
vertical resolution for spacebome radar can be well approximated by [45]: 


r eff * 'V <P h tan 4»i ) 2 + (y cos <f»j) 2 (3.5) 
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Essentially, equation (3.5) Imposes the same constraint on the antenna 
pointing angle as does equation (3.3). To achieve a small vertical resolution, the 
pointing angle of the antenna and the beam width also need to be small. With orbit 
heights of 300 km and 830 km, the swath widths are calculated In the following 
table 


Table 3.2 Swath Widths for Different Orbits and Antenna Pointing Angles 

Orbit Height swath width (<]> = 35°) swath width (<(> = 30°) 

300 km 420 km 346 km 

830 km 1162 km 923 km 


Because of the limits on the antenna pointing angle and the vertical resolution, the 
swath widths of the radar system are smaller than those of LAWS (9). 


CLUTTER FROM ANTENNA SIDE LOBES 

Clutter is generally defined as any unwanted radar echo. In the context of 
this dissertation. It mainly consists of the echoes from land and sea. The clutter 
from the main antenna lobe may not be a problem, since it can be separated from 
the received signal by range gating. However, the clutter from the antenna sldelobes 
can be a severe problem to the system. This is due to the fact that the backscatter 
coefficients of sea and land clutter are much larger than the volume -scatter 
coefficients of clouds (46-47). Even when the gains of the antenna side lobes are 40 
dB below the main lobe, the power of the clutter may still be higher than that of the 
signal returned. To solve this problem, we make some suggestions below: 

• Using a narrow pulse width and large antenna pointing angle, so that the 
energy of the clutter can be reduced as a result of reduced cross-section of 
clutter. 

• Designing the antenna such that sldelobes near vertical are very small. 

Although a detailed study of clutter rejection Is an Important issue In the RAWS 
study. It is beyond the scope of this dissertation. 
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PRF SELECTION 


With tp the pulse length, rj- the range resolution, and c the speed of light, the 
range resolution is r r = ctp/2. To obtain a fine range resolution, it is necessary to 
have a short pulse length. On the other hand, the frequency resolution requires the 
that minimum time of measurement be approximately equal to the inverse of the 
bandwidth of the Doppler filter. 

A typical radar pulse length is on the order of microseconds, and the target 
(wind) speed is on the order of tens of meters per second. To achieve a windspeed 
resolution of 1 ms' 1 , the required measurement-time is about 15 milliseconds at 3 
cm wavelength (10 GHz), and 4.3 milliseconds at 0.857 cm wavelength (35 GHz). 
Because of the required range resolution of 1 km, such long pulses cannot be used in 
the design of RAWS. Therefore, unlike a laser radar (lidar) system which may use a 
single pulse to measure the Doppler shift, a microwave radar wind sounder must use 
a train of pulses to measure the mean Doppler frequency and at the same time meet 
the requirements of range and frequency resolutions. 

To avoid range ambiguity for evenly spaced transmitted pulses, the pulse 
repetition time (PRT) is bounded by the following inequality: 

PRT > 2 T(j. ans + Techo + Tguard (3.6) 

where Ttrans 18 the pulse length of the transmitted pulse (it is expanded pulse width 
when the chirp technique is applied): T ec ho is the round trip delay of the 
electromagnetic wave propagating through a 20 km thick cloud; and T guarc j is the 
time required to compensate for variations of the Earth's surface. For a uniformly 
spaced pulse train, with Ttrans = 20 ps, T ec ho = 163 ps at 30 deg. antenna pointing 
angle, and Tg uarc j = lOps, equation (3.6) leads to the following inequalities: 


PRT > 213 ps 
PRF <S 4700 Hz 
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For a train of pulse-pair transmission as shown In Figure 3.4, and discussed later in 
chapter 5, the PRT Is bounded by 


PRT > 2T 2 + 2T trans + Techo + T gU ard (3.8) 

Assuming that the pulse width Is 20 ps, Ti=PRT, T 2 = 35 ps. and Tg uar( j = 10 ps. it 
follows that 

PRT > 283 ps 
PRF < 3523 Hz 

For a quadratic receiver, taking PRF as the sample frequency, the Nyqulst 
frequency is equal to half of the PRF. Thus, to measure 60 ms* 1 radial wind speeds, 
the minimum PRFs required to satisfy the Nyqulst criterion are equal to 

4v f 8000 Hz at f=10 GHz 

X = 128400 Hz at f=35.5 GHz 

From (3.9), we can see that there is a very severe frequency-ambiguity problem at 35 
GHz; the maximum Doppler shift can be as large as 8 times the Nyqulst frequency. 
Later in Chapter 5, we discuss the algorithms for reducing frequency ambiguity 
problems. 

3.2.2 ANTENNA SCAN SCHEMES 

Selection of the antenna scan mode was one of the major problems In the 
design of RAWS that we had to Investigate. During the measurements of Doppler 
frequencies, both the satellite and the antenna beams are in motion, and this may 
cause changes In the beam position and antenna pointing angles. Such changes 
may Introduce errors In estimation of the mean frequency as well as decorrelation 
to the returned signals. The changes In antenna pointing angle and beam position 
can be reduced or minimized by choosing an optimal antenna scan mode. In this 
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section, three scanning modes are examined: the first is uniform scan; the second is 
discrete scan; and the third is a combination of the uniform and the discrete scans. 
These three scan schemes require different levels of complexity in hardware and 
also introduce various amount of errors into frequency measurements. 

3.2.2.1 UNIFORM SCAN 

In a uniform scanning scheme, the antenna rotates with a constant angular 
speed achieved with a scanning motor or by satellite rotation. However, there is a 
severe problem with this simple scanning strategy caused by satellite motion and 
antenna beam rotation. This scan mode is illustrated in Fig 3.5a. For example, 
assuming the PRF = 3500 Hz and that each measurement contains 256 samples, it 
would take 256/3500=73.1 ms for each measurement. However, if the uniform 
scanning period T s =10 s, during one measurement interval the outer beam would 
travel 9.6 km and the inner beam 7.96 km with a 300 km orbit. These values are 
more than 10 times as large as one footprint of the antenna beams at f=35 GHz. 
Such a large displacement in antenna beam can cause a serious decorrelation in the 
returned signals. In addition, the antenna pointing angle changes nearly 2.6° 
during each measurement. Such a large change in antenna pointing angle is bound 
to cause a large error In calculations of the wind vectors. The only way to reduce the 
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Figure 3.5. Antenna scan pattern; a) uniform scanning; b) discrete scanning; c) a combination of the 
uniform and the discrete scanning. 


changes in pointing angle and displacement in antenna beam Is to lengthen the 
scan period. However, because of the required horizontal resolution of measured 
wind field, the density of measurements is limited to at least once per 100 km by 100 
km. This Imposes a limit on the scan period being less than 27 s. 

3. 2. 2. 2 DISCRETE SCAN 

The second mode, shown In Fig 3.5b, Is a discrete scan which assures that the 
pointing angle does not change during each measurement. The only decorrelation 
introduced Is due to the satellite motion; the beam travels 0.54 km during a 
measurement of 256 samples. However, there Is no antenna pointing error during 
the measurement. Strictly speaking, the discrete scheme needs to be semi-discrete, 
as the motion of the antenna beam from one position to the next cannot be at a very 
high speed. Otherwise, the phase lock loop (PLL) used to track the Doppler shift of 
clutter caused by satellite motion may lose track of the Doppler frequency from the 
clutter (see section 3.2.3). Considering the size of the antenna the discrete scan is 
difficult to Implement mechanically. Electronic scan for such a large antenna may 
also be too costly and complicated to be achievable. However, use of a set of phased 
arrays cannot be ruled out. 

3.2.2.3 COMBINATION OF UNIFORM AND DISCRETE SCAN 

The third scheme Is to focus the beam to the center of the resolution volume 
while a measurement is being made. This strategy can be Implemented by a 
combination of both electronic and mechanical methods. For example, the 
reflector of the antenna can be controlled by a motor, and the antenna feed can be 
an electronically scanned phased array. During each measurement, the feed moves 
the beam back and forth to focus on the center of the resolution volume, while the 
reflector rotates at a uniform speed. This scheme, shown In Fig 3.5c, can greatly 
reduce the decorrelation problem which was inherent In the uniform scan-scheme. 
However, the change In pointing angle, during each measurement, can introduce 
estimation errors. The angle change can be calculated by the following expression: 


where A<|> is the angle difference between the initial and final pointing angle in a 
measurement of the mean frequency, T m is the measurement time, R is the slant 
range, and v is the satellite speed. Assuming that the measurement time is T m =73 
ms (256 pulses), with satellite speed v = 7500 ms' 1 , then the angle difference Is equal 
to: 


{ 2 mrad for space shuttle orbit 
0.6 mrad for near polar orbit 

Such a small change in antenna pointing angles may not introduce any significant 
error in retrieving wind vectors as we will see later In section 3.3.3. 

In comparison of the three scan modes, we may conclude that the uniform 
scan is not a satisfactory scheme; the discrete mode may be too costly to implement; 
the combination of uniform and discrete scan is the optimal choice among the 
three. 

3.2.2.4 SCAN TRAJECTORY 

The coordinates of the two beam trajectories In a Cartesian system can be 
expressed as: 

x = h x tan 4>. x cos cot + vt 

y = h x tan x sin tot i = 1,2 (3. 1 1) 

where 

4>i = the antenna's elevation angle 
to = the angular speed of scan 
v = the speed of satellite 
h = the height of the satellite orbit. 

Figure 3.6 presents the trajectories of the antenna beams at orbits of 300 km and 
830 km. 
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Figure 3.6. Antenna scan trajectories at 30 deg. and 35 deg. pointing 
angles with scan period = 10 s; a) at 300 km orbit: b) at 830 km orbit. 








This figure shows that all four trajectories intersect only at very few points. It can 
be seen that if the measurements are made at the places where at least two 
trajectories intersect, the measurements would not be uniformly distributed. 
Namely the marginal area would get more measurements than the central region. 
The same problem has also been experienced in the research on LAWS. An adaptive 
laser shot pattern was suggested to cope with this problem. This adaptive shot 
pattern is designed to control slew rates of the scanner and schedule pulse 
suppression in regions of low Information potential 123] . A similar method can be 
applied to RAWS. In fact, the RAWS has the advantage of two antenna beams which 
can lead to a more evenly distributed measurement pattern. However, we leave this 
topic for future study. 

With the intention of simplifying the analysis of the error bounds of the 
estimation of wind speed in the next section, we will consider the measurements are 
being made only at places where two beam trajectories intersect. This may not be 
necessary in practice, since we have assumed that the wind field being measured is 
much larger than the resolution volume. Thus, all of the measurements of Doppler 
frequencies in a 100 km by 100 km area could be freely combined to derive the wind 
vectors in this area. However, measurements of Doppler frequency with pointing 
vectors close together may give more accurate estimates. 

To further simplify the error analysis In the following section, assume that 
we only want to make measurements at the intersections of the trajectories of the 
antenna beams with the same elevation angle. We can pre-calculate these points 
with the following equation 

(2n + l)x-2 0 = 0) t n = 0,1,2,... 

2 h tan<)>i cos 0 = vt 


or 


2coh tan<t>i 
v 


cos 0 + 20 -<2n + 1) re = 0 


n = 0,1,2,.... 


( 3 . 12 ) 


where 
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v = satellite speed 
h = orbit height 
<t>i = elevation angle 
0 = horizontal angle 

With T s =10 s, <f> 1 =30° and <}> 2 = 35°, the above equation can be solved numerically. 

Table 3.3 lists the positive angles at which two trajectories Intersect. Since the 
trajectories of beams are symmetrical about the ground track of the satellite, the 
negative angles can be obtained by adding minus signs to the values listed In Table 
3.3. 


Table 3.3 Roots of equation (3. 12) 


0 In rad ( <(> = 30°) 

01nrad(«|)=35 o ) 

0.307 

0.2567 

0.79 

0.702 

1.084 

0.961 

1.3357 

1.182 


1.38 


3.2.3 TRACKING THE DOPPLER SHIFT CAUSED BY SATELLITE MOTION 

In the context of calculating the error bound for wind speed we made the 
assumption that the Doppler shift due to the satellite motion could be accurately 
compensated. In the following, we will discuss a strategy used to compensate this 
Doppler shift. As shown In Figure 3.7, which Illustrates a simplified block diagram 
of the radar system, two stages are used In compensating the Doppler shift by 
tracking the mean Doppler frequency of clutter return. In the first stage, the 
compensation is accomplished through an inertial navigation system sending 
antenna pointing angle Information to a control processor. The control processor 
tunes a voltage-controlled oscillator (VCO) to counterbalance the Doppler shift 
caused by satellite motion according to the information given. This Is an open-loop 
system; the compensation error is dependent on the accuracy of the Information 
provided by the inertial navigation system. The second stage compensation is 
carried out by a closed loop tracking system. This tracking system can be 
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implemented by a digital phase lock loop, an analog phase-lock loop, or a 
combination of a Kalman filter and a phase lock loop. In the following, we assume 
that the second stage compensation is performed by a second-order analog phase 
lock loop. 



Figure 3.7. Functional Block Diagram of Radar Wind Sounder. It 
shows the two stage tracking of the Doppler shift caused by satellite 
motion. 


Denote the velocity of the satellite, u, the pointing vectors r and the Doppler 
shift fo. Asume the satellite is heading in x direction. Then u, r and can be 

expressed as: 
u = vx 

r= sin0 cos<t> x + sin9 sin«|> y + cos0 z 

and 

2u 

fD = "j[~ r (3.13) 
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Then, the bound of the remaining Doppler frequency after the first-stage 
compensation can be expressed as 

2 v 

Afp « -jj— [cos 0 cos 0 A0 + sin 0 sin 0 A0] 

2v 

< Y cos(0 - 0) max( A0. A0) (3.1 4) 

where A0 and A0 are remaining antenna pointing errors after the compensation by 
the first-stage tracking. If the pointing errors can be reduced to less than 1 mr by the 
inertial navigation system of the satellite (according to the LAWS report [91, the 
pointing error Is on the order of tens of microradians), the bound of Afp is about 2 

kHz at f=35.5 GHz. 

In the following, a computer simulation of the error of the dynamic tracking 
carried out by a second-order phase-lock loop is completed. Investigation of more 
advanced tracking methods is left to further studies. 


PD LPF 



Figure 3.8. Equivalent linearized baseband model of a PLL. 
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3.2. 3.1 FREQUENCY TRACKING BY A PHASE LOCK LOOP 


As shown In Figure 3.8, a phase lock loop (PLL) consists of a phase detector 
(PD), a low pass Alter, and a voltage -controlled oscillator(VCO). The low-pass filter 
can be first, second, or third order, which Is also referred to as the order for the PLL. 
The input signal to the PLL Is usually modeled as the sum of the signal s(t) and the 
noise n(t) 

y(t) = s(t) + n(t) 


with 


s(t) = A(t) sin((OQt + <t>(t)) 


n(t) = cos cogt + n s sin cogt 

where both nc and n s are narrow-band Gaussian noise processes. The phase 
detector can be modeled by an ideal multiplier whose output equals 

e(t) = K d [ sin(<D(t) + n'(t) )] 


with 


Kq - Kvco ^ 

K vco = amplitude of the VCO signal 
0 = phase of the input signal 
0 = phase of VCO 

A 

<t> = 0 - 0 = phase error 


and n'(t) is also a narrow band Gaussian noise which has the following relation with 
the Input noise n(t) 


n c (t) a n s (t) a 
n'(0 = — — cos0 +— t — sin0 
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LINEAR APPROACH 


When the phase error O is very small, the PLL can be approximated by a 
linear system with the following transfer function 


S?7=H<s> = 


„ „ F(s) 

Kvco Kq — - 


9(s) 


1 + Kvco Kq 


F(s) 


(3.15) 


For a second order PLL with perfect integrators, the error transfer function can be 
written as [48-491: 


e e (s> 

0j(i) 


2 2 
s +2£o) n s + co n 


(3.16) 


with ^ the loop damping factor, ©n the natural frequency. The input noise power to 
a phase lock loop is equal to 


+oo 

J S nn<“> d “- N 0 B IF 

oo 


(3.17) 


where S n n(w) Is power spectral density of input signal, Bif is IF bandwidth. The 
bandwidth of the phase-lock loop is generally much smaller than the IF bandwidth 
and is defined as: 

+oo 

B L=^ J I H nn (co) 1 2dco Hz (318) 


The SNR in the loop is defined as 
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(3.19) 


1 

P = — o = 


B 


2 

OA 

0 


n o b l 


= SNR; 


IF 


• B, 


where SNRi Is Input SNR For more detailed descriptions of these parameters in a 
phase-lock loop, one can refer to [48], [50], [51]. 

LOCK LIMIT. HOLD-IN RANGES 

The hold-in range is defined as the input frequency range over which the 
loop can sustain a lock status. The hold-in range for a frequency offset of PLL is 
shown [48] as 

Aco h = ± K V co K D (3.20) 


Similarly, the maximum permissible rate of change of input frequency is limited by 



(3.21) 


These two parameters are important in computer simulations of the PLL since the input 
signal cannot have a frequency offset or a change of frequency rate larger than those set by 
(3.20-3.21). Otherwise, the loop would lose the tracking. 


NON-LINEAR ANALYSIS 

Linear analysis is valid when the input phase error is small and the loop 
SNR is high. However, to study the dynamic behavior of PLL's, we need to consider 
the non-linear properties of the PLL. As an analytical solution is difficult to 
achieve for a second or a higher order PLL, we resort to computer simulation for the 
tracking problem by the PLL. First, we need to develop a state variable equation for 
a second order PLL. An equivalent function block diagram of a second order PLL 
with perfect integrators is shown in Figure 3.9. 
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Figure 3.9. Equivalent block diagram of a second order PLL with 
perfect integrators. 


The corresponding state variable equation is: 


d<t» / dt = - 2£co n [sin <& + n'] - x + d0/dt (3.22) 

2 

dx / dt = (d^ [sin <t> + n'J 


In the computer simulation of this stochastic system represented by (3.22), 
the following parameters are chosen: 

d0/dt = 2000 u(t) sinco s t, input signal to the PLL 
where u(t) is a unit step function. 

a)- = angular speed of scan 

c =0.5 
“n = 100 

n' = a narrow band Gaussian random process, generated by an 
AR(1) lllter. 
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Seconds 


b) 

Figure 3.10. Monte-Carlo simulation of tracking error by a second 
order PLL with perfect integrator; a) ensemble average; b) standard 
deviation of frequency error 
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A 4th-order Runge-Kutta method is used in the simulations. It was shown 
that using the Runge-Kutta method in computer simulation of a random system, the 
result is equivalent to the definition of Ito stochastic equation [52-54], Figure 3. 10 
presents the Monte Carlo simulation. As expected the RMS error of tracking 
frequency depends on the SNR in the loop. When the SNR is high the RMS error in 
frequency tracking is very small. 


3.3 MEASUREMENT ERRORS OF WIND SPEED DUE TO ANTENNA 
POINTING ERROR AND FREQUENCY MEASUREMENT ERROR 


3.3.1 BASIC EQUATIONS 


The Doppler frequencies observed by the radar can be represented as 


, 2 (u w- u s- u a> 
f Di= l r i 


i= 1,2,3 


where 


(3.23) 


u w = the velocity of the wind field 
u s = the velocity of the spacecraft 
u a = the velocity of the phase center of the antenna 
rj - the antenna pointing vectors 
X = the wavelength 

Assume that the Doppler shift due to the motion of the satellite and the rotation of 
the antenna can be accurately compensated, and let fi be the measured mean 
Doppler frequencies, be the antenna elevation pointing angles, 0i be the antenna 
azimuth pointing angles, and W x , Wy, W z be the x,y,z components of the wind 
velocity. We can rewrite equation (3.23) in a matrix form: 


F d = AW2A 


where 



f i 


w x 

II 

Q 

tbi 

f 2 

w = 

w y 


f 3 


w z 


(3.24) 


(3.25) 


rj = sin <j>icos 0j x + sin <)>j sin 0j y + cos <+)j z 
with i=l, 2.3 and J= 1,2,3 


(3.26) 
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So A has the following form: 


A = 


sin([)|COS0| 

sin(()2COS02 

sin^cos©^ 


sin<J>^sin0^ cos<j>j 
siiu^siitf^ cos^ 
sirut^sin©^ cos^ 


(3.27) 


In practice 4>i and 0j can only equal to a limited set of values which depend on 
the antenna pointing angles, the antenna scan schemes, and the measurement 
pattern used in the system. For example, the elevation angle fa can only be equal to 
30° or 35°. 


Furthermore, if we assume that the vertical component Is independent of the 
horizontal components and can be estimated separately from the horizontal 
components, we can write the above matrices in a two-dimensioned format. 


A = 


sin<(>jcos0i sin<j>jsin0| 
-sin<(>2cos02 sin^sii^ 


(3.28) 


and 


f d = 




W 

W = 

X 

f 2 

w y 


(3.29) 


Notice that the elements in F D of (3.29) are not direct measurements of the 

Doppler frequencies. They are equal to the subtraction of the Doppler shifts caused 
by the vertical component of the hydrometeor motion from the measured Doppler 
frequencies. 

3.3.2 ERROR ANALYSIS 

The major sources of errors in calculating the wind fields are a) estimation 
errors of the mean Doppler frequencies, b) estimation errors of antenna pointing 
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angles, and c) tracking errors of the mean Doppler frequency caused by satellite 
motion. In the following section, we assume that the Doppler shift caused by 
satellite motion is accurately compensated, and discuss only the error bounds of 
wind measurements resulting from antenna pointing errors and estimation errors 
of mean Doppler frequencies. 

3.3.3 ERROR BOUND CAUSED BY FREQUENCY UNCERTAINTY 

When we express the linear equation of (3.24) in a vector form, and under the 
assumption that there is no antenna pointing error, the error resulting from 
uncertainties in frequency measurements is equal to 

5W f = A' 1 5 F d X/2 (3.30) 

where 5F D is the frequency uncertainty and 6W f is the error in the wind-field 

estimates. To derive the upper bound of equation (3.30), we take A' 1 as a linear 
operator in a metric space. The definition of the norm of A' 1 in a metric space is 
equal to 

I|A _1 || =max || A^xH (3.31) 

ll*||= 1 

where x is an any unit vector in the metric space J55]. It is known for a linear 
operator in a metric space that the following relationship holds: 

II Axils || A || IMI (3.32) 

where x is any vector in the metric space. Using this equation, the upper bound of 
equation (3.30) can be expressed as: 

||8W f lU||A- 1 ||||5F D |U/2 (3.33) 


- 67 - 


Up to now we did not specify which norm or metric space to be used in this 
discussion. The choice of norm often depends upon the objective and analysis, as 
well as mathematical convenience. It Is usually difficult to find the norm of a 
linear operator in an arbitrary metric space. However, if the metric space Is chosen 
to be l^space, the norm of A' * is well derived as the maximum row sum of absolute 

magnitudes In A" *[55-56J: 



i = 1,2,...,M 


(3.34) 


where N and M are numbers of columns and rows in A respectively. With |[A~*|| 
defined as in (3.34), we can always find a unit vector 5 Fq such that when 5 Fd is 
substituted into equation (3.33), then the left side of (3.33) Is equal to the right side. 

Therefore, equation (3.33) is the minimum upper bound of measurement errors in 
the l^space. With values in a limited set of pointing angles 0 given in Table 3.3, 

from equation (3.34), we have numerically calculated the norm of A' 1 The results 
are listed in Table 3.4. 


Table 3.4 Norm of linear operator A' 1 in 1^ space 


01 in rad 


e 2 in rad 

JlAl'll 

0.307 

6.617 

0.2567 

6.866 

0.79 

2.843 

0.702 

2.69 

1.084 

4.27 

0.964 

3.06 

1.336 

8.57 

1.182 

4.6 



1.38 

9.2 


Therefore, under the assumed condition, from equation (3.33) the upper 
bound of the measurement errors of the wind fields due to uncertainty in frequency 
measurement is 
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|| 5W f || < 4.6 X || 6F d || 


(3.35) 


so that If the required wind measurement error Is to be less than 1 ms'*, the errors 
of mean Doppler frequency measurements need to be less than 23 Hz. This result 
may be overly optimistic since It does not take Into account errors In the 
measurements of the vertical components of wind vectors. 


3.3.3 ERROR DUE TO UNCERTAINTIES IN ANTENNA POINTING ANGLES 

The error analysis due to the antenna pointing error Involves the analysis of 
a non-linear operator In a metric space. The analysis of a non-linear operator Is 
usually much more difficult than that of a linear operator. If the operator is 
differentiable (in our case, it certainly is), the errors resulting from small errors In 
antenna pointing angle 0 and <|> can be expressed as 156): 

M N 

” 8W all-lH F D II < 33 « 

i-1 ' )-l ' 


where 5W fl are vectors of errors, and 50j are antenna pointing errors, and M and 

N are numbers of different elevation angles and different horizontal angles in 
matrix A respectively. Again, by using the relationship || A x || < || A || || x || and F D 

= A W 2/X, the above equation can be rewritten as 



where 
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(3.38) 


f M 

tA*'* 


V 


i=l 


N 



A 


50; 


J 


is a non-linear operator. We calculated the norm of this non-linear operator in a 1^ 

space with the angles listed in table 3.3. For 5<|> = 0.001 rad and 50 = 0.001 rad, the 
norm is listed in Table 3.5. 


Table 3.5 Norm of non-linear operator (3.38) in 1^ space 


Space shuttle orbit Near polar orbit 


01 In rad 

WHHHH 

01 in rad 

II- II 

0.307 

0.0086 

0.2567 

0.0093 

0.79 

0.0027 

0.702 

0.0029 

1.084 

0.0061 

0.964 

0.0034 

1.336 

0.0113 

1.182 

0.006 



1.38 

0.0125 


According to the value in Table 3.5, the wind error caused by an antenna pointing 
error of 0.001 rad can be calculated with equation 3.37: 

||5W a |U0.0125||w|| 

3.3.4 THE TOTAL ERROR BOUND 

The compound error in one measurement due to both uncertainties in 
frequency measurement and in pointing angle can be written as: 

||5W f ||sa||F D || + p||w|| (3.39) 

where a and p are constants whose values are dependent upon the antenna pointing 
angles, the scan mode, and the measurement pattern as well as the number of 
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Independent measurements. For a single measurement of a wind vector, under the 
condition discussed In the section 3.3.3, a is 4. 6 A. and p Is 0.0125. If the second term 
In equation (3.39) Is ignored, to obtain 1 ms' 1 accuracy In measurements of wind 
vectors, the frequency errors need to be less than 23 Hz. Notice, from equation 
(3.39), that the measurement error resulting from the antenna pointing error is 
increased with the true wind speed. 

3.4 CONCLUSIONS AND RECOMMENDATIONS FOR FURTHER STUDY 

a) The second term of equation (3.44) depends upon the true windspeed. When 
the windspeed increases, the error in the measurement due to a fixed antenna 
pointing-angle error also increases. 

b) Throughout the numerical calculations of the norms, we have assumed 
that the vertical components of the wind vectors were known. In practice, the 
estimation or measurement of vertical velocity certainly contains error. Therefore, 
the norm should be larger than those calculated in this chapter. 

c) We have assumed that the Doppler shift due to satellite movement can be 
very accurately compensated. Otherwise it will cause a bias in the wind 
measurement due to the compensation errors. 

d) The above error analysis is based upon the errors In measurements of the 
mean Doppler frequency and estimation of the Doppler shift of the satellite. 
However, these estimation errors and measurement errors are usually functions of 
the SNR of the received signal. A final analysis concerning SNR of the system 
should be performed. 

e) The error analysis conducted above is crude; there may be other factors 
which need to be considered in evaluation of the performance. A method used In the 
LAWS study may also be applicable to the RAWS: namely error analysis through 
computer simulation. In that method, the error is analyzed by using a Monte Carlo 
technique; each source of error Is considered independently, and values are assigned 
to the sources by using random number generators that have the appropriate 
distribution functions. This can be used as a topic of future study. 
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f) Comparing the prf of a Ildar system (about 8 Hz for LAWS) and the prf of a 
radar system (3480 Hz in our case), it seems that the radar system has an advantage 
of using higher prf than the Ildar system. However, for measurement of Doppler 
shift, the Ildar only needs one pulse while the radar system needs 64 or more pulses. 
Thus the radar system only has a slight advantage on prf over the Ildar system with 
the current technology. 
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Chapter 4 


Estimation of Moments of Power 
Spectrum 


4.1 MOMENTS OF SPECTRUM 

In Doppler radar signal processing, iare often required to estimates of the moments 
of the power spectra of radar echoes. An ith moment of the power spectral density of 
echoed signal, m^, is defined as: 


+oo 

Jf 1 S(f) df 

oo 

ni- =~ (4.1) 

1 +oo 

JS(f) df 

_oo 

where S(f) is the power spectral density function of the signal. For most weather 
radars, only the first three moments are important [441: 

(1) mean power of the signal or zeroth moment which indicates the total 
signal energy returned from the target. In wind measurements from clouds, the 
zeroth moment is proportional to the echo from hydrometeors, and it may also 
reflect the attenuation caused by clouds or precipitation. This information may be 
useful for deriving the precipitation rate. However, it does not directly relate to the 
winds. 


(2) the mean Doppler velocity or the first moment of the normalized power 
spectral density is a measure of the radial component of the target velocity. This is 
an essential parameter in the estimation of wind vectors, and reflects a weighted 
average radial velocity of wind fields. 
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(3) The second moment is useful for deriving the root-mean-square (RMS) 
width of the power spectrum (or spectral width). It is a measure of the velocity 
dispersion, i.e., the shear or turbulence within a resolution volume. 

In the study of airborne radar wind sensors, the first moment is the 
fundamental variable for deriving the Doppler shift caused by the wind-driven 
droplets in clouds. The spectral width derived from the first and the second moment 
may be useful for determining the error in estimates of the first moment. In this 
chapter, we discuss the algorithms for estimating the first moment of the power 
spectral density. Monte Carlo simulations were applied to evaluate the second- 
order statistics of the algorithms discussed in this chapter. Since all of the 
estimators are functions of the SNR and the spectral width, results of the computer 
simulations are presented with the RMS errors as functions of SNRs and spectral 
widths. The discussion and the computer simulations are based on normalized 
frequencies which range from -1 to 1. However, the results can be easily converted 
to a specific application like RAWS, by multiplying the results with the Nyquist 
frequency. 

4.2 INTRODUCTION TO RANDOM SIGNALS AND SPECTRAL ANALYSIS 
4.2.1 RANDOM PROCESS 

Before we discuss the algorithms for estimating the moments of the 
frequency spectrum, let us briefly review some pertinent concepts of random 
processes and their power spectral densities. 

STATIONARY PROCESS 

A process x(t) is said to be a wide sense(or weakly) stationary process if and 
only if both the expected value and the variance of the random process are constants 
and its covariance between the values at any two time points, ti and t2. depends 
only on (t 2 -t i ), the interval between the time points, and not on the location of the 
points along the time axis. Mathematically, we can express these conditions as: 
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E(x(t» = p 


var(x(t)) = ct 2 


R Xx (t 1 .t 2 )= ECxftj) x*(t 2 )) = R xx (t 2 -ti) 


(4.2) 


Radar echoes are oflen assumed to be stationary random processes; at least over 
short periods, most radar echoes do resemble stationary processes. There are cases 
In which the radar echoes cannot be considered as stationary processes. One such 
example is that the antenna sweeps over several beamwidths during a 
measurement, the radar echoes may come from several different resolution 
volumes. Therefore, the spectrum of the radar echo may be gradually changing 
during the measurement, and may not be considered as a stationary random 
process. 

ERGODIC PROPERTIES OF A RANDOM PROCESS 

In addition to the assumption of stationarity, we also assume that the radar 
echoes are ergodic random processes. A random process x(t) is said to be ergodic if 
its time averages equal its ensemble averages. In almost all of the cases, we may 
only be concerned with the ergodicities of a few parameters of a random process. 
For example, the ergodicities of the mean and the autocorrelation are essential to 
spectral analysis. 


4.2.2 POWER SPECTRAL DENSITY 

For a stationary random process, the autocorrelation function and the 
power spectral density are defined as 

R xx (t)= E(x(t+x) x *<*» (4.3) 


+oo 



oo 
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The autocorrelation function R^fx) Is the Fourier transform of the power spectral 
density S(f), and provides the basis for spectral analysis. 

As a practical matter, one does not usually know the statistical 
autocorrelation function. Thus, it Is often necessary to assume that the random 
process Is ergodic In Its first and second moments. The statistical autocorrelation 
function of an ergodic random process can be written as: 

+T 

1 i m 1 

R xx (x )= T-*» 2T / X < t+T > x *(‘> dt (4.5) 

-T 

From this equation, it can be shown that the PSD of an ergodic random process is 
equal to : 


S(f) = 


lim 

T-*o 



+T 

2' 

1 

2T 

Jx(t) exp(-j2xft) 

► 


-T 



(4.6) 


The expectation operator is required in the above equation since the ergodicity does 
not couple through the Fourier transform; that is, the limit in the above expression 
without the expectation operator does not converge in any statistical sense [57] (581. 


4.3 COMPUTER SIMULATION OF RANDOM SIGNALS 

4.3.1 SHAPE OF THE DOPPLER SPECTRUM 

There are two important parameters In describing a stationary random 
process; one is the probability density function; and the other is the power-spectral- 
density function. A radar echo consists of the sum of the return signal and noise. 
The noise is often assumed to be white and Gaussian. Since the radar echoes from a 
resolution volume are from multiple scatterers, it is reasonable to assume that the 
radar echoes have Gaussian probability distributions as implied by the central 
limit theorem. Hence, the power of a radar echo should have an exponential 
distribution. 
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The shapes of spectra of weather-radar echoes are Influenced by both the 
shapes of the transmitted pulses and the antenna patterns. If we assume that both 
the shapes of transmitted pulses and the radar pattern can be approximated by 
Gaussian functions, the shapes of spectra of weather echoes can be approximated by 
the convolutions of two Gaussian functions, and, should also resemble Gaussian 
functions. In addition, as the Fourier Transform of a Gaussian function results in a 
Gaussian function, the autocorrelation function of a weather radar echo is often 
assumed to be Gaussian. In the following discussion, we will always make the 
assumption that the shape of the spectrum of a weather radar echo is Gaussian. 

In the computer simulation of weather radar echoes, besides the shapes of 
spectra, we also need to consider the spectral widths of the radar echoes. In the case 
of RAWS, because of the high speed of the spacecraft, the spectral widths of the 
weather radar echoes are determined primarily by the beamwidth of the antenna. 
For example, as derived in Chapter 3, the 3 dB bandwidth of the Doppler frequency 
of a returned radar signal at 35 GHz is roughly equal to 

2v 2v 

B D=r p -r = 1875 H2 


where v is the velocity of the spacecraft, X is the wavelength, (J is the antenna 
beamwidth and L is the antenna aperture. The ratio of the bandwidth of the Doppler 
frequency to a 3500 Hz PRF is equal to 


B, 


1875 


2f N "3500 


= 053 


In a pulsed-Doppler radar, the PRF acts like a sampling frequency to the weather 
radar echoes. Therefore, from the above calculation, the bandwidth of the Doppler 
frequency of radar echoes is approximately half of the Nyquist interval. If we 
consider that turbulence, wind shear, motion of the antenna beam, and target drift 
all would add decorrelation effects to the returned signals and broaden their 
spectra, the actual radar-signal spectral widths may be greater than calculated 
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above. A quantitative study of the effect of antenna scan rate on signal 
decorrelation can be found In [59]. 


4.3.2 GENERATING RANDOM VARIABLES WITH SPECIFIC PDFS 

In computer simulations of radar signals. It Is often necessary to generate 
random processes with determined probability density functions and correlation 
functions. In this section, we review some of these methods for generating a random 
variable with a specified probability density function (PDF). In the next section, we 
discuss how to generate a random process with a specified autocorrelation function. 

Generating random variables with a specific PDF often starts with 
generating a series of uniformly distributed random numbers. Since the computer- 
generated random sequence is quasi-random, the random-number generator should 
be chosen according to the application. For Monte Carlo simulations with several 
hundred runs, and each run requiring thousands of samples, we need to choose a 
random number generator which produces statistically uncorrelated numbers with 
a very long cycle. For example. In the simulations performed in this chapter, the 
random number generator Is chosen from [60], and this generator has a cycle of 2^0. 

Generating random variables other than those with uniform distributions 
can be achieved through performing transformations on uniformly distributed 
random variables. It is known that If two random variables x and y are related by 
the expression x = fly), then the PDF of y is equal to [61): 


p(y) = p(x) 


dx 

dy 


(4.7) 


In particular, if x is a uniform random variable one can obtain the desired random 
variable by the transformation y=F" 1 (x), where F(y) is the distribution function of y. 
For generating some random variables, multi-variable transformations may be 
needed. Let xj and X 2 be uniform random variables and yj and y 2 be random 
variables with the desired probability distributions. Then the Joint PDF of yj and 
Y2 18 
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p(y\, y 2 ) dx dy = p(xj, x 2 ) 


g(X|, x 2 ) 

5(yi, y 2 ) 


dx dy 


(4.8) 


where 


3(xj, x 2 ) 
3(yi.Y2) 


is the Jacobian determinant. 


GENERATING A GAUSSIAN DISTRIBUTION 

An Important example of the use of variable transformations Is the Box- 
Muller method for generating random deviates with a normal distribution. 


p(y)dy = 


V2n' 


,-yW 


dy 


(4.9) 


Consider the transformation between two uniform random variables on (0,1), x^, 
and two quantities y lt y 2 > 


-2 In X| cos ( 271 x 2 ) 
y 2 = ^j-2 In Xj sin ( 2 XX 2 ) 


(4.10) 


Equivalently we can write 


1 ; y 2 

x 2 = -arctan- 


(4.11) 


Now the Jacobian determinant can readily be calculated: 

dfri. +y %)] 

^(yi. y 2 ) 


(4.12) 
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Therefore, and y^ have normal distributions as a result of the 
transformation. 

GENERATING AN EXPONENTIAL DISTRIBUTION 

The power of a weather radar signal is often exponentially distributed. It 
can be generated from a uniform distribution. Let x be a uniformly distributed 
random variable and y be an exponentially distributed random variable. With the 
transform x = F(y) = exp(-y). y has an exponential distribution. 

GENERATING A RAYLEIGH DISTRIBUTION 

The voltages of weather signals generally have Rayleigh distributions, 
which can be expressed as follows: 

P<y) = jJe'y 2/2(T ( 4 . 13 ) 

2 

Therefore, with the transform x = F(y) = e‘ y ^ 2a , y = V 2a(-ln(x)) has a Rayleigh 
distribution. 

4.3.3 GENERATION OF A RANDOM PROCESS WITH A SPECIFIC 
AUTOCORRELATION FUNCTION 

4.3.3.1 Using An ARMA Model 

A large number of random processes can be classified as certain 
autoregression and moving-average processes with p as the order of autoregression 
and q as the order of moving average. Such autoregression and moving average 
processes are usually denoted as ARMA(p,q). Even some random processes, which 
cannot be classified as ARMA(p,q), may be approximated by certain ARMA(p,q) 
models. Therefore, we may use an ARMA model to simulate the weather radar 
echoes. The spectra of radar signals returned from weather targets are expected to be 
Gaussian in shape. We may use an autoregression model. AR(p). to approximate 
such a process. 
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An ARMA(p,q) modeling only gives the desired autocorrelation function of 
the process. To make the random process have a specific probability distribution 
often requires a non-linear transformation after the ARMA modeling. In practice, 
such a non-linear transformation is often difficult to implement. A detailed study 
of such methods for generating random processes with specified spectra and 
probability distributions can be found in [62) . 

4.3.3 .2 Using An Inverse Fourier Transform 

Here, we describe a method using an Inverse FFT to generate weather-like 
signals. This method was discussed in [63] and [64J. We applied this method In 
generating weather radar echoes with Gaussian power spectral densities and 
exponential probability distribution functions in the computer simulations of this 
chapter. As discussed in section 4.2.2, a weather radar signal normally has a 
Gaussian-shaped spectrum, and its power is exponentially distributed. The detailed 
steps for generating weather signals with an inverse FFT method are depicted below: 

• The first step is to specify the power spectrum of the weather radar echo 
from a target by 

2 2 

G„ = -p= — e~^ n " fm> /2a n=0,l/.«,M-l (4.14) 

n V2jca 

where 

fn = n th discrete frequency 
M = total number of discrete frequencies 
Gn = the discrete spectral coefficient corresponding to fn 
f m = the power-weighted mean frequency 
a = is the standard deviation of the spectrum, which is defined as 
spectral width. 

In all the computer simulations of this chapter, f n and fm were normalized with a 
Nyquist interval [-1.0, 1.0). When f n or f m are not in the range of -1.0 to 1.0. aliasing 
occurs and the spectrum folds into the interval [-1.0, 1.01. 
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• Second, we assume the additive noise present In the real signal is white, i.e., 
it has a constant spectral density, and also consider that both simulated signal and 
noise powers must be exponentially distributed to represent the weather signal. 
Then, we can write the power spectral coefficient for frequency f n as: 


Sn= -ln(x n ) 


KGn+-^i 

M 


(4.15) 


where 0 5 x n 51 and is uniformly distributed. Pn. the total noise power, can be 
arbitrarily set to unity. K Is equal to 



(4.16) 


• The third step is to decompose the power spectrum into its real, A n , and 
imaginary, B n , components, adhering to the requirement that the phase angles of 
the sinusoid comprising the spectrum (i.e., the phase spectrum of the time signal) be 
uncorrelated and uniformly distributed between (-ji.ji). 

* - 1/2 „ 

A n = S n cos(2icy n ) 

1/2 

Bn = S n sin(2rey n ) (4.17) 

where y n is a random variable uniformly distributed between(0,l). Notice that S 1 / 2 
has a Rayleigh distribution, while A n and Bn have Gaussian distributions. 

• The final step is to generate the complex time signal by performing an 
inverse Fourier transform of A n and Bn. 


M-l 

z n = I n + jQ n =X (A i + i B i )ei2,tni/M (4.18) 

i=0 

where M is the number of samples. 
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Figure 4.1. Example of computer simulated spectrum of a random 
signal with spectrum width a = 0.3, mean frequency -0.7, and SNR 10 
dB. 

An example of computer simulated weather radar echo is shown on Figure 
4.1. The frequency is normalized to [-1,11. The mean frequency and the standard 
deviation are -0.7 and 0.3 respectively. In the computer simulations of this chapter, 
the total number of sample, M, is chosen as 1024. 

4.4 ESTIMATION OF THE MOMENTS OF THE DOPPLER SPECTRUM 

The algorithms for estimating of the first moment of the power spectral 
densities of radar signals can be classified into the following following three 
categories: 

• Algorithms based on the Fast Fourier Transform(FFT) of the signal. 

• Algorithms based on the covariance function of the signal (also called 

pulse-pair method). 

• Algorithms based on parametric modeling, such as the autoregressive (AR) 

and the moving average (MA) models. Only the autoregressive method 
will be discussed in this chapter. 

The algorithms based on the FFT method and the ARMA-model method 
require evaluation of the power spectral density (PSD) of the return signal first, and 
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subsequently the mean of the PSD. The algorithms based on covariance function do 
not require the calculation of the PSD of the signal. The ARMA methods may not be 
as efficient as the covariance and the FFT methods when the order of the model Is 
high. In the following sections, we will discuss each of these algorithms in more 
detail. 

4,4.1 FFT METHOD 

The amplitude spectrum for a digital signal is often expressed as: 

oo 

X<e)“) = Ix(n)e'i“ n (4.19) 

"OO 


In practice, however, the amount of data is always limited. One of the most widely 
used PSD estimators, based upon an FFT operation. Is typically referred to as the 
periodogram. For data samples xo,...,xn-i. the perlodogram estimate of the PSD is 
defined as: 


A 1 
S(f)- 


NAt 


N-l 

At X x n e ' i2nfnAt 


n=0 


2 


(4.20) 


where N is the number of data. At is the sampling Interval, and f is in the range of - 
l/(2At) < f < l/(2At). Use of the periodogram permits us to evaluate the PSD at N 

equally spaced frequencies f m = mAf Hz, for m=0,l N-l and Af= 1/NAt. If the Af 

factor is incorporated into S(f), then (4.20) can be written as: 


S =S(f )Af = ^ IX I 
m m N I ml 


2 


N-l 


e -j2jcmn/N 


n=0 


N 


(4.21) 


where X m are the coefficients discrete Fourier transform (DFT) of x(n). 
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Not only Is the perlodogram a biased estimator of the power spectrum, but it 

A 

is also not a consistent estimator, i.e. when N goes to infinity, S(f) may not converge 
to its mean value statistically [58]. This phenomenon is primarily caused by the 
absence of an expectation operator in the above equation. Several methods can be 
used to reduce the variance of the perlodogram PSD estimators. One of these 
methods is to divide a long data sequence into M short sequences, and separately 
applying each sequence to the perlodogram. The results are then averaged. 


Estimation of the mean of the PSD with the perlodogram is straightforward. 


M-l 


'f S. 


m m 


1 m=0 


f F M-l 


X S m 

m=0 


(4.22) 


1 i a 

where fp = is the Nyquist frequency, and defined as - fp. f ranges from 

-1.0 to 1.0. It can be shown that (4.22) generally is not a consistent estimator for the 
mean frequency. Only when the spectral lines are mutually independent, is (4.22) a 
consistent estimator (see Appendix 4.1). 

In equation (4.22), we did not consider the effect of noise. However, the 
estimate of mean frequency of the PSD Is often biased by the noise. If we write the 
power spectral density for Sffm) as a sum of the spectra of the signal and noise, 

S(fm) = + n m 


with Sjn and n m as the discrete spectral densities of the signal and noise 
respectively at frequency fm. Equation (4.22) can now be expressed as follows: 
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M-l 

Yf(S m +n m ) 
^ m m m 

* 1 m=Q 

f p M-l 

y (S m + n ) 
m m 

m=0 

M-l 

y/m (S m +n m ) 

1 m=0 

"f F (S + N) 


m s m N 
= 1 + N/S + 1 +S/N 


(4.23) 


where m s is the first moment of the signal spectrum, m N is the first moment of the 

noise spectrum, and N/S is the inverse of the SNR Therefore the bias caused by the 
noise is equal to 


m s N/S mjvj 
B,as = 1 + N/S + 1 + S/N 


(4.24) 


Notice that in equation (4.24) the bias is a linear function of m s When the noise is 
white, the second term of (4.24) is zero, andthe bias approaches zero as the mean 
frequency, m s> approaches zero. This is illustrated by Figure 4.2. When the mean 
frequency is zero, there is no bias caused by noise. 

In addition to the bias caused by noise, the estimator is also biased by 
aliasing of the spectrum of the signal. From Figure 4.2, it can be observed that, at 0 
dB SNR, the estimate suffers a considerable amount of bias. The estimate is about 
half of the true value of the mean frequency. At 10 dB SNR, the bias Is negligibly 
small in the middle part of the Nyquist interval, and the estimates almost agree 
with the true values of the mean frequencies. However, when the mean frequency 
approaches the ends of the Nyquist interval, the estimates starts to curve towards 
the 0 frequency. This phenomenon is caused by the aliasing of the spectra of the 
signals. In summary, the FFT estimator of (4.22) only works well when the 
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-1.0 -0.5 0.0 0.5 1.0 


Mean Frequency (Hz) 

Figure 4.2. Simulation of FFT estimator without de-allaslng and 
noise suppression, under the condition a = 0. 1, SNR = 0 dB, and 10 dB, 
and the number of runs = 200. 

signal-to-noise ratio Is high and the signal is free of frequency aliasing. In case of 
poor SNR and frequency aliasing, the estimated mean frequencies yield very large 
errors as shown in Figure 4.2. In the next two sections, we discuss noise suppression 
and de-aliasing algorithms that can remove or reduce the bias due to the noise and 
frequency aliasing. 

4.4.1 .1 Estimation of The Mean With Noise Suppression 

If we assume that the signal and noise are uncorrelated, and the spectral 
density of the noise can be estimated separately, one way to reduce the bias caused 
by frequency aliasing is to subtract the noise spectral density, N(f m ), from the 
derived spectral density and calculate the mean of the resulting spectrum, i.e., 

M-l 

? 1 m=0 / a «p\ 

w - 25) 

m=0 
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-1.0 -0.5 0.0 0.5 1.0 


Frequency (Hz) 
b) 

Figure 4.3. Monte Carlo simulation of FFT estimator with noise 
suppression; a) estimate under the condition a = 0. 1. SNR = 0 dB. and 
10 dB, with 200 runs; b) RMS errors as functions of frequency and the 
spectral width of power spectrum of signal. 
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Although this estimator may result in some spectral lines with negative values, the 
biases are effectively removed from the estimates. This method was examined by 
Sirmans and Bumgarner [64] with computer simulations. Their results can be 
summarized as follows: the FFT with noise suppression is an estimator of mean 
frequency unbiased by noise even for a low SNR and the standard deviation of the 
estimate with the noise suppressed is comparable to the standard deviation of the 
signal plus noise mean estimates. 

Some of the results of the computer simulations for the noise suppression 
scheme are presented in Figure 4.3., which shows that the bias caused by noise is 
removed by the noise suppression. However, the estimates of the mean frequency of 
the signal are still biased at the frequencies near the ends of the Nyquist interval 
because of frequency aliasing. 

4.4.1. 2 FFT Method with De-allasing 

As shown in previous sections, the performance of the FFT method is 
severely degraded by aliasing of the spectrum. This Is particularly noticeable when 
the spectral width of the signal Is large. The effective unambiguous frequency 
measured by the FFT estimator can be significantly reduced without application of 
a de-aliasing method. One method to reduce the effect of aliasing is to shift the peak 
of the spectrum to the center of the Nyquist interval (at zero frequency) before 
applying the FFT estimator, and add the shifted frequency back in the final 
estimate. Such an algorithm was discussed by Zmic [65). In the computer 
simulations performed In this chapter, we used a similar method as explained 
below: 

• a) Smooth the periodogram by weighted running average. The size of the 
window used in the smoothing was selected proportional to the spectral 
width of the PSD. The maximum size of the window was limited to half of 
the Nyquist interval. 

• b) Search for the peak of the smoothed spectrum of the data, then locate Its 
position, say fp. 
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• c) Shift the entire spectrum such that the peak fp Is now at zero frequency. 

• d) Apply the following formula to estimate the mean frequency , 


M-l 


I 


f S' 
m m 


f - f p + 


1 m=0 
f F M-l 

S S m 

m=0 


where S'm - Sm-p- 


(4.26) 


In Implementing such an algorithm, tradeoffs are encountered In choosing 
the size of the smoothing window. If the window size Is too small compared with the 
spectral width of the data, the algorithm may select a frequency peak which Is not 
close enough to the true mean frequency, and the bias caused by frequency aliasing 
cannot be effectively removed. On the other hand. If the width Is too large, it will 
increase the uncertainty of the location of the peak frequency. In the computer 
simulations performed, we chose the widths of the smoothing windows to be 
proportional to the spectral width of the signals. 


Some results of Monte-Carlo computer simulations of such an algorithm are 
presented In Figures 4.4. In these computer simulations, the Nyquist intervals are 
normalized to [-1, 1]. The standard deviation of the spectrum are chosen to be 0.1. 
The RMS error of the estimate Is defined as: 


RMS ERROR 



(4.27) 


A 

where f m is mth estimate of the mean frequency f, and M is number of runs. 
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- 1.0 - 0.5 0.0 0.5 1.0 


Mean Frequency (Hz) 
b) 

Figure 4.4. FFT Estimator with de-aliasing applied with a=0.1, 
SNR=10 dB, a) the plot of ensemble average of estimates; b) plot of 
RMS error with 200 runs. For the RAWS parameters, the scale would 
be from -1.75 kHz to + 1.75 kHz. The frequency shown here is 
normalized to the Nyquist frequency. 
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These results show that the de-aliasing algorithm improves the performance 
of the FFT estimator significantly. Even at a SNR of 0 dB. the estimate is still able 
to converge to the true mean frequency, indicating the de-allased FFT estimator is 
an unbiased estimator of mean frequency. 

The fact that the de-aliasing method reduces the bias caused by noise is in 
agreement with equation (4.25). This equation indicates that, as the mean 
frequency goes to zero, the bias caused by noise also goes to zero. After applying step 
c) in the de-aliasing algorithm, the bias caused by noise acts like a zero-mean 
random variable. Hence, the estimates can converge to the true values of the mean 
frequencies. However, this is only true when the noise is white For colored noise, 
noise suppression may need to be applied before applying the de-aliasing algorithm. 

4.4.2 THE COVARIANCE ESTIMATOR 

The covariance method for computing the moments of the Doppler spectrum 
has come into widespread use in recent years. These methods have been discussed in 
many papers: Rummler[66], [67), Benham and Groginsky [68], Miller and 
Rochwarger[69), Sirmans and Bumgamer[64[, [70J, and others. The covariance 
method is a time-domain estimator. Therefore, it does not need to estimate the 
power-spectral density. One obvious advantage of this method is that it requires 
fewer computations. In addition, this method does not require equal time interval 
between sampled pairs. This property makes it possible to combine this method 
with waveform modulation for removal of range and frequency ambiguities (see 
Chapter 5). 

The following covariance method for estimating the mean of the spectral 
density is described by Sirmans and Bumgarner [64J. This method, also called 
pulse-pair processing, is based on the fact that the moments of a random variable x 

mn= E(x n ) (4.28) 

are related to the derivatives of its characteristic function 4>(a)). the Fourier 
transform of its probability density function with a reversal in sign, by 
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(4.29) 


•n 


m n ~ 


d n 0(0) 
d x n 


Since the autocorrelation function and the power spectral density S(f) constitute a 
Fourier transform pair, the moments of the power spectral density are related to the 
derivatives of the autocorrelation function by an equation similar to (4.29): 


, n R'(0) 
m n = ’ R(0) 


(4.30) 


Expressing the covariance function in a polar form 
R(x) = A(t) exp[j27ig(x)] 

where A(x) is a real even function of x and g(x) is a real odd function of x. The mean of 
the PSD is equal to 

H[H] -M 

x=0 x=0 


Notice that A'(x) is an odd function so that A'(O) =0. This is also true for g(0). 
Therefore f can be approximately written as 

sW^s<0> _ SM _ i Arg[R(t>] (4.32) 

T TKX 

In the computer simulations, a maximum-likelihood unbiased estimator of R(x) was 
used [69]: 

a N-l 

R* x (x) = l/N £x n+ iX* n W.33) 

n=0 
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(4.34) 


A 

f = 


1 

— arctan 
nx 


f Im(R 


^Re(R 


A 

XX 

~ 

XX 


(X)) 


(X)) 


\ 

J 


where X n+ i and Xn are complex samples of data, spaced x seconds apart, and f 
ranges from -1.0 to 1.0. 

One Important fact to note Is that a large Interspacing between two pulse 
pairs does not reduce the accuracy of the covariance method. On the contrary, 
according to [71). when the Interspacing increases, the accuracy of the covariance 
method should also Increase because more independent samples of R(x) can be 
achieved. The results of the Monte-Carlo simulation presented In Figure 4.5 agree 
with this fact: the estimator with longer spacing between pulse pairs give smaller 
RMS errors. The data used were computer simulated weather radar signals with 
Gausian spectra. To select different spacing of pulse pairs, the following formula 
was used to calculate Rfx): 

M 

ft(x) = X x(i*k + 1) x*(i*k) 

1=1 

where k is interspacing of pulse pairs. Figure 4.5a shows that the covariance 
estimator is a consistent estimator of the mean frequency. 
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Estimate (Hz) 




Figure 4.5. Computer Simulation of the Covariance Method with 128 
pairs of samples a) mean estimate, with lnter-pulse-pair spacing = 1, 
7, SNR = 10 dB. a = 0.1; b) RMS errors with 200 runs. For the RAWS 
parameters the scale would be from -1.75 kHz to + 1.75 kHz. The 
frequency shown here is normalized to the Nyqulst frequency. 
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4.4.3 THE PARAMETRIC METHODS 


Many deterministic and stochastic discrete-time processes encountered in 
practice are well approximated by certain rational transfer function models. Such 
models, known as the ARMA models, represent random processes, whose PSD are to 
be estimated in terms of linear difference equations of the following form: 

P 3 

x(n) = X 3p,i x(n- i) + ^ bq ^ e(n-k) (4.35) 

i=1 k=0 


where bqo is 1 and e(n) is a zero-mean white Gaussian noise sequence. x(n) may be 
viewed as the response of a linear time-invariant filter whose input is a white noise 
sequence. The transfer function of such a filter has the form 


H(f) 


X^ j2 * fk 

6(f) k=0 


4>(f) 






k=l 


(4.36) 


Equation (4.36) consists of two parts: the autoregressive part and the moving- 
average part. The autoregressive (AR) portion consists of the poles of the filter and 
is the denominator of H(f); the moving average (MA) portion consists of the zeros of 
the filter and is the numerator of H(f). The PSD of the random process x(n) is given 
by 


S xx (0 = H(-f) H(f) S nn (f) = i H(f) I 



(4.37) 


where <* nn is the variance of the input noise, or the power spectral density of the noise. 
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Figure 4.6. Computer Simulation of AR(2) and AR(20) estimators 
with 200 runs, SNR = 10 dB. cr=0. l f and 128 samples for each run; a) 
plot of mean estimates; b) plot of RMS errors. For the RAWS 
parameters the scale would be from -1.75 kHz to + 1.75 kHz. The 
frequency shown here is normalized to the Nyquist frequency. 
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To examine the performance of the parametric method in evaluating the first 
moment of the power spectrum of radar echoes, we will concentrate our attention on 
the Autoregressive models. Because the power spectrum of the radar echo is expected 
to be Gaussian in shape, it can be approximated by an all-pole form of a transfer 
function, or an Autoregressive model as follows: 

P 

x(n) = I a p/i x(n-i) + e n ( 4 . 38 ) 

i=l 

where 

e n = white noise 

p = number of poles of the autoregressive model 
ap.i * coefficients 

There are two ways to solve for the coefficients a p i: one Is to use the Yule-Walker 
equation: the other is to use the Levinson-Durbin Algorithm. Both of these methods 
are discussed in numerous papers and text books concerning time-series analysis 
[57), [72], [73] and [74]. In the following, we will only discuss the Levinson-Durbin 
algorithm. The Levinson-Durbin algorithm requires only order O(p^) operations, 
as opposed to 0(p 3 ) for Gaussian elimination in solving the Yule-Walker equation. 
The algorithm proceeds recursively to compute the parameter sets : 

2 2 2 
{ an, CT-j), {a2i, a 2 2, o 2 ), ..., {a pl , a p 2, ..., a pp , cr p } 

where aij are estimates of the coefficients, and a p is the estimate of the variance of 
e n - Note that an additional subscript has been added to the AR coefficients to denote 
the order. The final set at order p is the desired solution. In particular, the recursive 
algorithm is Initialized by 
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RxxO) 

311 " Rxx(O) 


®1 = 0- I an | 2 )R xx (0) 

with the recursion for k = 2, 3 p given by 


a kk = - 


k-1 

+ X a k- 1, n^xx(k' n ) 

; Q=i 


’k- 1 


(4.39) 


a ki = a k- 1. i + a kk = a k- 1, k- 1 


(4.40) 


CT k = (! *l a klJ 2 ) CT k- 1 (4 41) 

Once the coefficients have been calculated, the power spectral density can be 
determined from the following 


S A r(D = 


q 2 At 

P 

1 + £ 
k=l 


(4.42) 


The estimate of the first moment can be determined from equation (4.22). Some of 
the computer simulation results are presented in Figure 4.6. Figure 4.6a shows the 
ensemble average of the estimates for AR(2) and AR(20): Figure 4.6b shows the RMS 
errors for AR(2) and AR(20). It shows that the increase of the order of the AR model 
does not improve the accuracy of the estimates. 


4.4.4 RANDOM SAMPLES AND SPECTRUM ESTIMATION 

It is well known that with the conventional equally- spaced samples, 
aliasing will occur if the power spectrum falls out of the Nyquist interval. However, 
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with unequally spaced samples, we may achieve estimates of the power spectrum of 
a stationary random process that are not aliased. Perhaps the best known result Is 
that Poisson sampling of a stationary process is alias-free [75]. In the following 
section, we choose a sampling scheme 

t n = nT + y(AT n - 03)T (4.43) 

where AT n Is a random variable having a uniform distribution on Interval [0.1] and 

y is in the range [0,1]. To evaluate the power-spectral density, we used an Intuitive 
formula, which Is similar to the definition of PSD for a continuous random signal, 
given below 


N-l 

S(0 = X x(t i ) e ’ j(0ti A i (4.44) 

n=0 

where A^ = t j+ j - 1,. To evaluate the first moment, we used the definition: 

1/2B 

m =B™oo! J fS<0df «•«> 

-1/2B 

When y = 0, this sampling scheme is identical to the equally spaced sampling 
scheme and the power-spectrum estimate is equivalent to the DFT method. When y 
= 1.0, this sampling scheme achieves maximum randomness. Computer 
simulations were performed using this scheme. The results show that for a signal 
with a small spectrum width (on < 0.1), the random sampling scheme does give an 
alias-free estimate of the power spectrum with reduced SNR. However, when the 
spectrum width is large, the shape of the estimated PSD is not discemable. Some 
work done on random sampling can also be found in [76], 
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4.5 RESULTS OF COMPUTER SIMULATIONS, CONCLUSIONS, AND 
FUTURE WORK 

4.5.1 MONTE CARLO SIMULATION OF THE ESTIMATORS 

To compare the statistical properties of the algorithms discussed in this 
chapter, several Monte-Carlo simulations were completed, each with 200 runs and 
128 samples, to calculate the RMS errors for each of the algorithms under different 
frequencies, spectral widths, and SNRs. The frequencies are normalized to two 
Nyquist intervals; one is the unit Nyquist Interval from -1.0 Hz to 1.0 Hz, the other 
is from -1750.0 Hz to 1750.0 Hz used by the RAWS. The RMS errors are calculated 
according to equation (4.27). 

The results are presented in Figures 4.7 to 4. 12. In Figures 4.7 to 4.9, the RMS 
errors are plotted as functions of both SNRs and frequencies. In Figures 4. 10 to 4. 12, 
the RMS errors are presented as functions of SNRs and spectral widths of the 
signals. The results Indicate that all three estimators are unbiased; the ensemble 
averages of the estimates converge to the true mean frequencies of the signals. The 
RMS errors of estimates in these algorithms are independent of frequency. 
However, the RMS errors are linear functions of the spectral widths of the spectra of 
the signals In the range of 0.0 Hz to 0.4 Hz for the unit Nyquist interval, or in the 
range of 0.0 Hz to 700 Hz for RAWS’s Nyquist interval. 

In comparing of Figures 4.7 to 4.9, it can be seen that when the SNR is high, 
for example 10 dB, all three estimators give similar performance in term of mean- 
square errors in estimation. However, when the SNR decreases to 0 dB, the 
covariance method produces the smallest estimate errors of the three; at 0 dB SNR, 
the AR(2) estimator is slightly better than the FFT estimator. Previous work 
performed with these estimators appeared to indicate that the AR(p) modeling 
should require more computation than the FFT estimator. However, when the order 
of the autoregressive model p is low, such as 2 in our simulation, the estimator 
based on an autoregressive model may need less computation power than the FFT 
estimator. 
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- 1.0 - 0.5 0.0 0.5 1.0 

-1750 -875 0 875 1750 

Mean Frequency (Hz) 
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- 1.0 - 0.5 0.0 0.5 1.0 

-1750 -875 0 875 1750 

Mean Frequency (Hz) 
b) 

Figure 4.7. Monte Carlo Simulation of FFT estimator with de- 
allasing. with 200 runs, and 128 samples for each run; a) Expectation 
simulated; b) RMS error as a function of SNR, and frequency. Both 
normalized (-1, 1) and RAWS (-1750, 1750) frequency scales are 
shown. 
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Figure 4.8. Monte Carlo Simulation of the Covariance estimator, the 
spectral width a =0.3, with 200 runs, and 128 samples for each run; a) 
Mean estimate for SNR = 0 dB; b) RMS errors as functions of SNR and 
frequency. 
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b) 

Figure 4.9 Monte Carlo Simulation of the Autoregressive, AR{2), 
method, with 200 runs, and 128 samples for each run; a) Expectation 
simulated; b) RMS error as a function of SNR and frequency. 


- 104 - 




Spectral Width (Hz) 


Figure 4. 10. Monte Carlo simulation of FFT estimator with RMS 
error as functions of SNR and spectral width. A fixed mean 
frequency 0.5 is used In this simulation. 



Spectral Width (Hz) 


Figure 4.11. Monte Carlo simulation of the covariance estimator 
with RMS error as functions of SNR and spectral width. A fixed mean 
frequency 0.5 Is used In this simulation. 
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Spectral Width (Hz) 


Figure 4.12. Monte Carlo simulation of AR(2) estimator with RMS 
error as functions of SNR and spectral width. A fixed mean 
frequency 0.5 is used in this simulation. 


4.5.3 CONCLUSIONS AND FURTHER WORK 

The de-aliased FFT, the covariance estimator, and the AR estimator are 
comparable in performance. The errors in the estimates depend on the spectral 
width and SNRs of the input signals. The errors are mainly functions of the spectral 
widths of the signals when the SNR is above 5 dB. As discussed in the introduction, 
the spectral widths of the signals in RAWS are expected to be half of the Nyquist 
interval. In such cases, all the estimators discussed in this chapter produce large 
estimation errors, about 10% of the Nyquist frequency. This means, for a 3500 Hz 
PRF, the RMS errors are around 175 Hz. Part of the estimate error is caused by 
aliasing of frequency spectrum. To reduce the estimation error caused by aliasing of 
the spectrum, we need either to increase the PRF or to reduce the inter- spacing 
between the two pulses in the covariance method. In Chapter 5, we discuss how to 
use the latter combined with waveform modulations to mitigate the ambiguity 
problems as well as to reduce the estimation errors caused by aliasing of spectrum. 
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In this chapter, we discussed the performance of several estimators of mean 
frequency: the FFT estimator, the covariance estimator, and the estimators based 

on autoregressive models. The covariance estimator seemed to produce slightly 
smaller errors than the FFT estimator and the autoregressive-model-based 
estimators in the computer simulations. However, this is true only under the 
condition that the power spectrum is symmetric and has only one peak. In practice, 
Interference caused by leakage from the transmitter, or clutter, may cause the 
returned signals to be non-symmetric or have more than one peak. In these cases, 
the covariance method may not work as well as the other estimators. We only 
discussed two autoregressive-model-based estimators, AR(2) and AR(20). The results 
showed that there was no difference between these two estimators In terms of the 
estimation errors. We also discussed the random sampling method to estimate 
spectra of signals. Because of limits on time and volume, we did not investigate the 
random sampling method In great depth. 
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APPENDIX 4A PROOF THAT THE FFT ESTIMATOR FOR MEAN IS 
CONSISTENT 

To prove that (4.22) is a consistent estimator of mean frequency, we must 
show that as M approaches infinity, the variance of (4.22) approaches zero. Let S(fj) 
be denoted as Si Af, and without loss of generality, assume that 


M- 1 

Xsw = i , 

i=0 
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(4.46) 


Theoretically, the above result is bound by the following inequalities: 


M-l 2 M-l M-l M-l M-l 

X ffCTj Af 2 I X Vj P s s CT.CT.Af 2 < X X fif| CT.CT.Af 2 = (E(mi )) 2 

1=0 1=0 J=0 J » J J i=0 J=0 J J 


where 


(4.47) 


mi = the first moment of the power spectrum 
M = total number of samples 

PLJ = cross correlation coefficients between spectral line i and spectral 
line j 

CTj = standard deviation of spectral line i 
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When the spectral lines are mutually independent, the cross-correlation 
coefficients are 


Pi,j = 


P 

lo 


when i = j 
when i * j 


Under this condition, when M goes to infinity, the variance of the FFT estimator 
approaches zero because Af *» 1/M. Therefore, the estimator is a consistent 
estimator. On the other hand, if the spectral lines are not mutually independent, 
the variance is bounded by the square of mi- Therefore, when spectral lines are not 
mutually independent, only when m^ is equal to zero is the FFT estimator a 
consistent estimator. 
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Chapter 5 


Algorithms for Removal of Range and 
Frequency Ambiguities 


5.1 INTRODUCTION 

The characteristics of radar echoes from weather targets impose limitations 
and tradeoffs on applications of Doppler radar systems for weather observation. 
Such limitations come from two facts: a) weather targets are distributed quasi- 
continuously over large spatial regions (tens to hundreds of kilometers), and the 
strength of radar echoes from a significant weather target easily spans an 80 dB 
power range [59]; b) the mean Doppler frequencies of radar echoes from weather 
targets are often higher than the maximum unambiguous frequency of the radar 
system. In other words, the inherent ambiguity problem in radar systems becomes 
more prominent in weather radars. 

For pulsed-Doppler radars, the unambiguous range, r a> Is generally defined 
as the maximum distance that a transmitted pulse can travel to a target and echo 
back to the radar receiver before the next pulse Is transmitted. The maximum 
unambiguous frequency is generally considered as half of the PRF. When a target 
area is located beyond the unambiguous range, the echoes returning from that area 
arrive after the next pulse is transmitted. This creates what is commonly referred 
to as range ambiguity. Radars having uniform PRF and without some form of 
coding from pulse to pulse cannot discriminate between echoes coming from targets 
located within the unambiguous range.and those outside this range. 

The maximum unambiguous velocity, v a . measured by a Doppler radar is 
related to the maximum unambiguous frequency by 


v a = ± X prf 14 
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where X is the radar wavelength. A velocity higher than the maximum 
unambiguous velocity causes frequency ambiguity. Range ambiguity and velocity 
ambiguity are not independent; it can be shown that v a r a =Xc/8. This relation 
indicates that either the unambiguous range or unambiguous velocity can be 
increased only at the expense of the other. 

For example, in the conceptual model of RAWS described in Chapter 3. the 
PRF is set to 3500 Hz to avoid range ambiguity. However, if the radar needs to 
measure 60 ra s" 1 wind speed, the maximum Doppler frequency of radar echoes 
caused by wind will be 14200 Hz at 35 GHz. This frequency Is 8 times as large as the 
Nyquist frequency of 1775 Hz. and will certainly create frequency-ambiguity 
problems. In such a case, the true mean frequency of the signal cannot be measured 
correctly with the conventional methods because of frequency aliasing. On the 
other hand, if we choose the PRF high enough to satisfy the Nyquist criterion, say 
PRF = 29000 Hz, there will be no frequency ambiguity. However, with a maximum 
range of 20 km the radar echoes from different transmitted pulses will overlap. 
Therefore, to satisfy the required maximum frequency and range of RAWS we have 
to solve the radar ambiguity problems. 

The algorithms developed thus far for reducing radar ambiguities can be 
classified into one of the following categories: 

• Interpulse phase coding: in an interpulse phase coding method, the 
transmitted pulses are modulated with a sequence of discrete phase codes. At the 
receiver, a coherent reference signal Is used to correlate with the received signals 
from a specific trip. The objective of the phase coding is to make signals returned 
from different trips have poor cross correlations. Therefore, the interference 
between echoes from different trips can be reduced. In general, after coherent 
processing, the interference from overlaid echoes appears like white noise. 

• Multiple pulse-repetition-frequency methods or multiple pulse-repetition- 
time methods. These methods are also known as staggered PRF or staggered PRT 
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methods. In these methods, pulses are transmitted with two or more different PRF's 
or PRTs. The estimates from each PRF or PRT can be combined or correlated to find 
the true mean frequencies of the signals. With staggered PRF or staggered PRT 
methods, the maximum unambiguous frequency of a radar system can be extended. 

• Polarization coding from pulse to pulse and frequency hopping from pulse 
to pulse: Doviak and Sirmans [81] suggested an orthogonal polarization coding for 
successive pulses. In this method, two orthogonal polarizations are used for two 
successive pulses. Therefore, the overlaying between the radar echoes from the first 
trip and the second trip is reduced by the polarization of the antenna. However, the 
depolarization by hydrometeors and the radar system limits this method to about 
20 dB of suppression of the interference. The frequency-hopping method was 
suggested by Doviak and Zmic [44]. In that method, consecutive pulse pairs are 
transmitted at different frequency steps: therefore, the echoes from different pulse 
pairs can be separated by filters. 

• Since the wind profiles are normally continuous in both frequency and 
space, the frequency ambiguities could be corrected by removing the discontinuities 
in the wind profiles. Such a method was demonstrated by Jiro [82], 

The first three algorithms mentioned above are related to inter-pulse coding 
in either phase, polarization or pulse position. In this chapter, we will further 
discuss these methods, and compare their performance with computer simulations. 
Since the radar ambiguity function is a widely used tool for analyzing and studying 
the ambiguities of waveforms of radar signals [83], we will also discuss the 
algorithms for removal of radar ambiguities in term of their radar ambiguity 
functions. To compare the performances of different algorithms, Monte Carlo 
simulations were performed to calculate the second-order statistics of these 
algorithms. In addition, a new algorithm for reducing range and Doppler 
ambiguities using waveform modulation is discussed. 
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5.2 INTERPULSE PHASE CODING 

In the following sections, we will examine the performance of different 
ambiguity-removal methods through the Interpretations of their radar ambiguity 
functions. As we know, a radar ambiguity function can be used to show the 
properties of range resolution, frequency resolution, and the distribution of 
ambiguities for a particular waveform. We will correlate the radar ambiguity 
functions of the algorithms discussed in this chapter with their performance In 
estimating the mean frequencies. Before discussing the algorithms, we briefly 
review the properties of ambiguity functions. A detailed discussion on ambiguity 
functions of radar waveforms can be found In [82]. 

5.2.1 RADAR AMBIGUITY FUNCTIONS 

The response function for a radar signal has two basic forms: 


+oo 

X^T, 4>) = J u(t) u*(t+T) e‘i 2n4>l dt (5.1) 

-OO 


+oo 

X^T, <(>)=/ U(f+«t»)U*(0 ei 2nfx df (5.2) 

-oo 


where u(t) is the transmitted signal and U(f) is the Fourier transform of u(t). The 
response functions given in (5.1) and (5.2) were derived by using the definition of a 
matched filter and thus are sometimes called the matched-filter response functions. 
The ambiqulty function is defined as I % 1 2 . 

The function xft.^) for a fixed 4» and a fixed x describes the amplitude 
modulation of the signal at the output of a receiver filter from a target with a 
Doppler shift <j> relative to the center frequency of the filter and a delay x from the 
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time to which the filter is matched. The response function can also be interpreted as 
the correlation between the transmitted complex waveform u(t) shifted by the 
Doppler frequency u(t)e‘J TO t >t , and u(t) itself, where the shift <(>=0 occurs at time To- 
When <(»=0, the response function reduces to the autocorrelation function of the 
transmitted signal u(t). To avoid misinterpretation, hereafter in this chapter we 
refer to I x I as the ambiguity function. 

Several of the important properties of a radar signal can be determined from 
its ambiguity function. An ambiguity function has its peak value centered at the 
origin, indicating that the largest signal output occurs when a target has the range 
and velocity to which the filter is matched. In practice, being matched to a 
particular range and Doppler means that the filter Is (1) sampled at the time 
corresponding to the round-trip delay of the transmitted signal to a target and (2) 
tuned In frequency to the Doppler shift corresponding to the radial velocity of the 
target. 


Targets which appear at ranges and velocities such that I x(t. 4>) I Is about as 
large as I x(0,0) I are indistinguishable to the radar. The width of the peak about (0.0) 
defines the resolution of the waveform. Other peaks away from (0,0) correspond to 
the ambiguities of the waveform. The minimum distance in the delay-time domain 
between the peaks other than the origin and the peak in the origin corresponds to 
the maximum unambiguous range; similarly, the minimum distance in the 
Doppler-frequency domain between the peaks other than the origin and the one in 
the origin corresponds to the maximum unambiguous frequency. 

Ambiguity functions can also be used to study clutter rejection. Clutter can 
be described as any unwanted backscatter. Waveforms having x(t,4>) =0 in the region 
of the x<J> plane where clutter exists generally have good clutter-rejection properties. 
More precisely, the summation over the entire x<> plane of the product of I x(t» 4>) 1 2 
with the clutter distribution over the plane determines the total interfering 
clutter signal. The energy of the clutter can be calculated from the following 
equation: 
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(5.3) 


C = JJ p(x, <» o(x, <{>) I^( x * <W| dx d0 

-oo 

where ct(x, 4>) is the backscatter cross section of the clutter, and p(x,(t>) is the density 
function of clutter. Notice the spreading loss factor (4 tc)3 r 4 is omitted in the above 
equation. 

A cross-ambiguity function describes the situationwhen the receiver filter Is 
matched to a modulation v(t) which is different from the transmitted modulation 
u(t). The cross-ambiguity function is defined in the following forms: 


+oo 

X uv (x, * J u(t) v * (t + x)e*^ 7t< ^ )T dt (5.4) 

.OO 
+ OO 

or X uv (t. <t>) = J U(f+<|>> V*(0 e' 2nh df (5.5) 

_oo 

A cross-ambiguity function can also indicate information about the range and 
frequency ambiguities as does an ambiguity function. 

As examples of ambiguity functions, the following figures present two 
ambiguity functions: one is for a rectangular pulse: 

Itls i 

' ® lO elsewhere 


where 5 is the pulse width. The ambiguity function of u(t) is 
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The ambiguity function of above u(t) is 


f /x\ jrccbx ~ \ 1 1) sin + 4>K5 — I x I ) 
_ ] rect \25/ e 5 jc(<xt + 4»)(5 -Ixl) 


0 when I x I > 8 


As we have mentioned, the ambiguity function represents the output of a matched 
filter. The Doppler frequency information is contained in the factor e^ T . These 
ambiguity functions are plotted in Figure 5.1. As indicated by equation (5.7), the 
peak of a signal with a Doppler shift would also shift In time by amount xo=<t>o/a. 
However, as a is usually very large, on the order of io + *^> to wou ^ ver y sma ^ 111 
comparison to 5 and can be ignored in most cases. 


5.2.2 RANDOM-PHASE CODING 

Random-phase coding is a method used for removing range ambiguities. A 
random-phase method was first applied in the RONSARD radar [84] using the 
random phase inherent in the magnetron transmitter. The phase coding was used 
to remove the bias in the estimate of the mean frequency due to overlying echoes by 
measuring the transmitted phases and using them to compensate for the phases of 
the received signals. At the receiver, only the echoes from the first trip were made 
coherent, and all the overlaid echoes from other trips were incoherent and appeared 
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as Increased white noise. Therefore, the bias in the estimate of mean frequency- 
caused by Interference could be reduced. 

Laird suggested a phase-coding scheme In which the phase of the transmitted 
pulse Is coded with binary random number 0 or re [85]. By employing a coherent 
reference, this scheme could re-cohere radar echoes from any particular trip and 
make the echoes from other trips appear like white noise. However, because the 
phase coding method only spreads the energy of Interference echoes in the 
frequency domain, and does not reduce the total energy of the Interference echoes, 
only one trip can provide a sufficiently large SNR for reliable measurement. As a 
result, the random phase coding method has limited use. 

Therefore, the random phase-coding scheme cannot effectively recover the 
echoes from a trip if the energy of interfering echoes from other trips is strong. For 
this reason. Siggla developed an adaptive filter for processing the random phase 
coded radar signal [86]. The adaptive filter can reduce the effect of interfering 
echoes and increase the effective SNR However, a significant improvement can be 
achieved only when the spectrum of the radar echo has a narrow width [87], and this 
method is computationally intense. 

Computer simulations of ambiguity functions for three different phase- 
coded pulse trains were presented in Figures 5.2, 5.3 and 5.4. Figure 5.2 shows the 
ambiguity function of a uniformly spaced pulse train with no phase coding. As 
expected, the uncoded pulse train generates periodic range and frequency 
ambiguities. Figure 5.3 shows the ambiguity function of a uniformly spaced pulse 
train with random phase coding as used in [85). We can see that the range 
ambiguities are reduced, the energy of interferences is distributed across the whole 
Nyquist interval and acts like white noise. However, the frequency ambiguities are 
unchanged. A similar statement can be applied to Figure 5.4. which shows the 
ambiguity function of a uniformly spaced pulse train modulated with a 13-bit long 
Barker-code. 
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Figure 5.1. Examples of ambiguity function: a) ambiguity function 
for single rectangle pulse, b) ambiguity function for single linear FM 
modulated pulse. Absolute value shown for response function. 
Ambiguity function is square of response function, but contrasts are 
too great to show effectively on a 3-D graph. All other graphs in 
Chap. 5 labelled “ambiguity function" are actually response-function 
magnitudes. 
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Figure 5.3. Ambiguity function of uniformly spaced pulses with 
random phase coding, the total pulse number is 30. 
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Figure 5.4. Ambiguity function of uniformly spaced pulses with 
Barker phase coding, the total pulse number is 13. 


In summary, random phase-coding methods can be used to correlate signals 
returned from a specific trip and make the echoes from other trips appear like white 
noise. Thus, the frequency measured will not be greatly biased by the interference 
from other trips. Taking another point of view, the phase coding methods can be 
viewed as multiplying the transmitted pulses by a random signal sequence in the 
time domain; therefore, without a matched receiver, the spectrum of the received 
signal would look like white noise too. 

However, the random phase coding does not reduce the total energy of the 
interference. Instead, it spreads out the spectrum of the interference over the entire 
Nyquist interval. Therefore, the interferences would act like white noise and 
deteriorate the SNR if proper filtering were not applied. Since only the signal 
returned from one particular trip would have large enough SNR ( > 3 dB) for reliable 
estimate of the mean Doppler frequency, in practice the random phase-coding 
scheme may not be able to retrieve the mean frequency from echoes of an arbitrary 
trip. 
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5.2.3 DETERMINISTIC PHASE CODING 


Zmic and Mahapatra further studied the adaptive filter method discussed by 
Siggia for processing random phase-coded signals, and concluded that effective 
Improvements in suppression of overlaid echoes are possible only when overlaid 
echoes have narrow spectral widths [87], Sachidananda and Zmic proposed an 
alternative phase coding method which can reduce the correlation of an overlying 
echo signal to zero at one-pulse lag [88]. With this phase coding, the covariance 
("pulse-pair") estimator can give an unbiased estimate of mean frequency in 
presence of overlaid echoes from adjacent trips. The sequence of codes they 
suggested is n/4, -7c/4, x/4.... At the receiver a sequence of -jc/ 4, 0, -x/4,... is used to 

correlate the first trip echoes, and a sequence of 0, -re/4,0 is used to correlate the 

second trip echoes. This phase-coding scheme can make the autocorrelation 
functions of either the first trip echoes or the second trip echoes equal to zero. 
Therefore, the measurements of the mean frequency for either the first pulse or the 
second pulse are not affected by echo overlaying, provided that the overlaying of 
echoes is only due to two adjacent pulses. 

The ambiguity function of a pulse train with such a coding scheme is shown 
on Figure 5.5. Figure 5.2 shows the ambiguity function of a corresponding pulse 
train without phase coding. One interesting point to note here is that one can arrive 
at a conclusion similar to that drawn from [88]. That is. all of the peaks of the range 
ambiguities from even trips are reduced to zero along the x (4>=0) axis. This is 
expected since the ambiguity function is reduced to the autocorrelation function 
along the x axis. However, we also found that this scheme shifted the peaks of the 
ambiguities along the x =±T, ± 3T... axes towards the x axis for all of the even trips. 
This means that the unambiguous frequency for the second trip echoes are reduced 
to a smaller range. 

From the above discussion, we can conclude that all the phase coding 
methods can be used to reduce the range ambiguity which results In flattening the 
peaks of the ambiguity function on the x<|> plane. However, the distances between 
these peaks along the frequency axis is not changed by the phase coding methods. 
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As a result, the maximum unambiguous frequency remained unchanged. Therefore, 
phase coding is suitable for high PRF situations where the frequency ambiguity is 
not a problem. 



Figure 5.5. Ambiguity function for uniformly spaced pulses with 

-jc/ 4, rc/4 -ic/4, n/4 inter-pulse coding.the pulse shape is Gaussian. 

and total number of pulses is 30. 


One other problem with the phase-coding methods is that the leakage from 
the transmitter may interfere with the receiver's operation if the isolation between 
the transmitter and receiver is not nearly perfect. This may be a serious problem 
for a chirp radar as an expanded pulse often lasts several tens of micro-seconds, and 
during the transmitting period the receiver cannot receive useful signals. 


5.3 MULTIPLE PRF AND FREQUENCY AMBIGUITY REMOVAL 

There have been many staggered PRF or PRT methods proposed for solving 
range and frequency ambiguities. Hennlngton (1981) suggested a method of 
transmitting a train of short pulses following by a long pulse to reduce the 
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frequency ambiguity of the targets 189). The long pulse was used to estimate the 
reflectivity of the targets; the short pulses were used to estimate the moments of the 
frequency spectrum. However, this scheme may not be applicable to a spread target. 
Sirmans et al suggested a staggered-prf method using two or more prfs to extend 
unambiguous frequencies (90). Similar methods also were proposed by 
Sachidananda and Zmic [881 and Ludloff and Minker (91). In the remainder of this 
section, we will examine only two of these methods: the method discussed by Ludloff 
and Minker will be called "STAGGERED PRF METHOD-A." and the method 
discussed by Sachidananda and Zmic will be called "STAGGERED PRF 
METHOD-B." 

5.3.1 STAGGERED PRF METHOD — A 

This algorithm was intended to be used in solving the blind speed problem of 
a MTD ( moving target detector) radar [91], but it may also be applied to Doppler 
weather radar applications. The algorithm resolves the frequency ambiguity by 
illuminating a target with a set of pulse bursts with different PRF's Fk. where 

A 

k=l,2 K. Let fo,k he the measured Doppler frequencies associated with prf Fk. 

These measured frequencies may be different from the true frequencies because of 
frequency aliasing. However, for each prf Fk the true frequency can be found among 
the following fj k frequencies: 

= *b,k + J k Fk i=— / -1,0,1,... (5.8) 


To search for the true frequency, the algorithm systematically searches for those 
integers Ik that cause all estimates fj t k with different k to fall within a small 
frequency interval or correlation bin. The average of the estimates provides an 

A 

improved estimate fp with reduced standard deviation 
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(5.9) 


, K 

%=k £ ( Ww 

k=l 


and extended unambiguous frequency because fr may be much larger than any Ffc. 
From now on, let us consider the case that K = 2 with the stagger ratio 


F]:F 2 = m:n 


where m and n are relatively prime numbers and m < n by definition. The expanded 
unambiguous frequency interval is 

f u = nF 1= mF 2 (5.10) 


PROBABILITY OF FALSE CORRELATION 

Although this algorithm reduces the standard deviations of the estimates of 
the mean frequency, it does introduce false correlation occasionally. The false 
correlations occur mainly at those locations where fi.k frequencies approach one 
another with the minimum possible distance dmin hi frequency. There is a simple 
relationship between dmin and Ffc given by [91]: 


d — — F 
d min “ n = 2 F 1 
n 


(5.11) 


Let us assume the estimates are normally distributed. Then, if we define the 
difference between the two estimates as x, the probability density function for x has 
the following form : 
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expf- 


(5.12) 



2 2 
I 2 (a! + ct 2 ) 

V J 


where xls equal to 0 for a true frequency measurement and equal to d m j n for false 
measurements. Since 



if we let F 2 = 1, the false measurement probability can be approximated as: 

05/n 

P((f - f'p)> 0.5 durin) =1 - J p(x) dx 

-0.5/n 


0.5/n 


= 1 - 


1 

f ^ ) 

, exp 

OjV 2x(l + (m/n) 2 ) 

J 

2o?((m/n) 2 + 1) 

V J 


dx 


(5.13) 


-0.5/n 


Here, <ri is approximately a linear function of the spectral width of the signal. As a 
result, when the standard deviation of the power spectral density increases, the 
probability of false measurement also increases as shown in Figure 5.6. It also can 
be observed that the difference between the false measurement and the true 
frequency is usually one or more Nyquist intervals. Thus, it is possible to correct 
these false estimates by continuity or the method we are going to discuss in the next 
section. The fi.k can be estimated with the FFT method, the AR method, or the 
covariance method. 
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b) 

Figure 5.6. Simulation of STAGGERED PRF-A method, a) a = 0. 1, SNR 
= 0 dB, b) a = 0.3, SNR = 0 dB. Normalized frequencies are shown. 
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5.3.2 STAGGERED PRF METHOD - B 


This staggered PRT method was described by Zmic and Mahapatra [87]. In 
this scheme, the pulses are transmitted alternatively with two different pulse 
repetition intervals, Tj and T 2 . The frequency estimate is determined by a formula 
derived from the covariance method: 


f 1 f R(T 2) >i 

r T~ ji(T 2 - Tj) ar ° ^R(Ti) J 

1 /-R(T2)x 

"tiT 2 (1 -K) arS ^R(T!) J 


(5.14) 


where R(t) is the autocorrelation function of the radar echoes, and K is a constant 
equal to T 1 /T 2 . The maximum unambiguous velocity depends on the difference of 
Tj and T 2 , and is equal to: 


f 

‘max ~ . T|) 


(5.15) 


The variance of this estimator can be derived from (5. 14), 


VAR(f) = 


jt(T 2 - T t ) VAR(0 1 ‘ V 


(5.16) 


where 0i and 02 are equal to arg(R(Ti)) and arg (R(T 2 )). A more detailed analysis can 
be found in [87]. From this equation, we can observe that the standard deviation of 
the estimate increases as the difference between Ti and T 2 decreases. The results of 
a computer simulation are presented in Figure 5.7. These results indicate that the A 
method is superior to the B method in terms of RMS errors: however the B method 
does not give false estimates (individually big errors). 
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Mean Frequency (Hz) 


Figure 5.7. Simulation of STAGGERED PRF-B method with T 1 /T 2 = 
0.875, cr=0.3, SNR = 0 dB. Normalized frequency Is shown. 


5.3.3 RANDOM PULSE POSITION CODING 

The multiple PRT methods discussed above can be classified as special cases 
of random pulse position coding. The general pulse train with pulse-position coding 
and phase coding can be expressed as: 

N-l j<b (t-nT-A) 

y(t) = £u(t-nT 0 -A n )e n n (5.17) 

n=0 

where A n is the displacement from nTo and <|>n represents the phase coding. The 
form of the ambiguity function for equation (5.17) and the properties of the 
ambiguity function of the position coding were reported in [921. We will not proceed 
further on this topic. 
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5.4 WAVEFORM MODULATION AND AMBIGUITY REMOVAL 

5.4.1 WAVEFORM MODULATION REDUCES RANGE AMBIGUITY 

The current algorithms developed for removal of range and frequency 
ambiguities are primarily Involved with inter-pulse coding in either phase or time. 
In this section, we will discuss the utilization of waveform modulation on 
transmitted pulses to reduce range ambiguities. Furthermore, as shown later, 
waveform modulation can be combined with the covariance method discussed in 
Chapter 4 to solve the frequency-ambiguity problem. 

Waveform design has found many applications in radar and communication 
systems. In radar systems, waveform design is often related to clutter rejection and 
pulse compression. In communication systems, it is used In the area of spread 
spectrum communications for code multiplexing [93 - 94]. 

The basic concept of waveform modulation for range-ambiguity reduction 
can be illustrated in an example from code multiplexing In a communication 
system. Waveform modulation allows two different signals to be transmitted in the 
same bandwidth and time interval, and at the receiver, the two signals can be 
separated by the process of correlation detection. This can be shown 
mathematically: assume s i (t) and S2(t) are two different waveformsfor codes), and 
sift) + S 2 (t) is the receiver input. For receiving signal sift), a reference signal sift) is 
used at the receiver to correlate with the input signal. Let T be the signal integration 
time, and the correlator's output may be expressed as: 

T 

y(t) = ^ J (s 1 (0 + SjfOJsjft) dt 
0 

T T 

1 f 2 1 

=Y I sjft) dt + y J SjftJSgft) dt (5.18) 

0 0 
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The first term In (5.18) Is the signal output and the second term represents 
Interference from S 2 (t). The Interference level Is therefore determined by the cross 
correlation of si(t) and S 2 (t). The usefulness of such a system Is determined by the 
cross correlation of the codes used. There are various coding techniques, such as 
Barker codes, m-sequence, and random phase-coding. Most of these codes are able 
to keep the sldelobe at the level of 1/N of the main lobe, with N Is the total length of a 
code (94). 

Instead of using discrete coding (mostly binary), we may consider using 
continuous waveform modulations. The main task is to find a set of continuous 
waveforms with low cross correlations and study the effective SNR at the receiver. 
One such candidate Is linear frequency modulation (FM). also known as chirp. For 
example, two pulses modulated with different llnear-FM waveforms may have low 
cross correlation. To simplify the following discussion, we consider waveforms 
having the same time Intervals and bandwidths. There may be many other linear 
or non-linear FM waveforms with low cross correlations. However, the discussion 
of those waveforms is beyond the scope of this study. 

Assume two Gaussian pulses u(t) and v(t) are modulated with conjugate linear 
FM, namely the first pulse In the pair is modulated with up chirp and the second 
with down chirp: 

u(t) = (2a) 1 /4 e' rtat2 + iJtal2 (5.19) 

v(t) = (2a) 1 /4 e" rtat2 ' )7tctt2 (5.20) 

where (2a) 1 / 4 is a constant to normalize the maximum value of the ambiguity 
function to 1. The second assumption made Is that the pulses are transmitted In 
pairs spaced by a time Interval T2. The transmitted pulse-pair train Is expressed as 
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(5.21). 


N-l 

y(t) = ^ua-nT^H- v(t-nT r T 2 ) 
n=0 

At the receiver, there are two channels. One is matched to the signal u(t) and the 
other is matched to v(t) as shown in Figure 5.8. The output from channel I and 
channel II are given as correlations between the received signals and the impulse 
response of the matched filter. 



Figure 5.8. Configuration of dual-modulation receiver. 

To examine the output of each matched filter, we need to calculate the 
ambiguity function and cross-ambiguity functions of the waveforms represented by 
(5.21). In fact, the cross-ambiguity functions are the outputs of the matched filters. 
The ambiguity function of such a pulse train can be calculated from the following 
integral: 
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+oo 


x u (x,<)>) = Jy(t)y (t + x) dt 


= X Su (*,<t>) + X Lv (x.«t>) e ' i27t<t>T 2 + x £uv (x - T 2 , 4 >) 
+ Xivu (x + T 2’ < f )e ' J27l<l>T2 


The analytical forms of x £u . X £uv .X £v and X £vu are listed below, 
derivation of these functions is included in Appendix 5.A. 


N-l 


X 


sin((N-m)TQ7t<t>) 

e F(N-l-m)4.T 0 x u (t-mT 0 ,4,) 


‘0« 


* Eu I when | x-mT 0 | < 5 


^0 


elsewhere 


r N-l 


L 


sin((N-m)TQ7C(j)) 

-j7 C (N-l-m)<)>T 0 x v (x-mT 0> <|» 


X S v (t ^) 


m=-(N-l) 

when | x-mTQ | < 8 


'-0 


elsewhere 


(5.22) 

A detailed 
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r n-i 


z 


e ^-'-“'♦ T 0 * vu ( , -mV)- 5 J S r nt) 


sin((N-m)TQ 7 t 4 ») 


W)= m sW ) 

Ivu I when | x-mTjj | < 5 


vQ 


elsewhere 


r N-i 


z 


e -j Jt (N-l- m )<t,T oZuv (T-m V ) 


=i 


m=-(N-l) 

when | x-mTg | 5 5 


^0 


elsewhere 


sin((N-m)TQTt<t») 

5 « 


Also from appendix 5 A. we obtain 


. , .J 2 ( ax +<|>)‘ 1 
X u (x.*)= + — T“ 


( 5 . 23 ) 


X v (x.« = e -” ,2 ( at2 * l 2 T l!L ) 


(5.24) 




2 9 

7t , 2 a<l> .re , .2 _ . a 4 > , 

e *2 ' aX + 2 2 e'*2 ‘“ l ' 2< * >T + 2 2 (5.25) 

a + a a + a 
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X uv (t,4>) 


W 


a + jot 


n , 2 

e'2 ' aT + 


^^(-a.^S+x 


a 2 + a 


a4> 


_2 . „2 
a + a 


( 5 . 26 ) 


The outputs from channel I and II are the cross correlations: 

channel I: Xe u M>) + Xj; U v^ -T 2’^ ^ 5 ' 27 ^ 

channel II: X£ V M) + %£ VU (* + T 2’^ (5.28) 


In equations (5.27) and (5.28), the first terms represent the signal outputs 
from the matched filters and the second terms represent interference. The 
ambiguity functions and cross-ambiguity functions given In equation (5.27) and 
(5.28) are plotted in 3D drawings in Figures 5.9 and 5.10. It can be observed that the 
energy levels of the Interference at each channel are much lower than that of the 
main lobe. In fact, it can be shown that the Interference levels are approximately 
equal to the inverse of the time-bandwidth (BT1 products of the waveforms. 

For a single target, the effective ratio of signal and interference can be 
calculated by the following equation: 



/a 2 + a 2 

Xe U v(°’°) 

\ ^ 


( 5 . 29 ) 


We are going to prove that the above expression is equal to the time bandwidth 
product of a pulse. Let the chirped pulse width be 8. The bandwidth then is 
approximately equal to a8, and the time bandwidth product, BT, is equal to 

BT = 8 a5 
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The parameter a and the pulse width 5 are related as follows 


1 



If we substitute these equations into (5.29), the effective signal to Interference ratio 
becomes 



( 5 . 30 ) 


Equation (5.30) above shows that the effective SNR for a single target is 
determined by the time bandwidth product of the FM modulated signal. By choosing 
large time bandwidth products for the transmitted pulses, we can reduce the 
interference level to an arbitrarily small value. 



Figure 5.9. Ambiguity function of dual linear FM modulated pulse 
train. T 2 = 0.3 Tj, time bandwidth product = 10.0 
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b) 

Figure 5. 10. Cross ambiguity functions, a) as an output of channel I, 
b) as an output of channel n. 
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However, In weather radar the target is always spread. The above 
calculation of effective SNR cannot be used for a spread target. The interference can 
be calculated from equation (5.3). However, from observing the plot in Figure 5.10, 
the domain of integration can be simplified to the area represented by the product of 
the time duration (T2-6 ,T2+5) and frequency Interval (-oo t +oo) since echoes only 
from this area can contribute to the interference. Denote C as the voltage output of 
the total Interference, it can be shown that (Appendix 5-A): 

+OO 

J X uv ( i_ T 2^) dxd<J» 


(5.31) 

Therefore, the effective signal to interference ratio is equal to 
s Ix(O.O)! 2 V a 2 + ot 2 BT 

C 2" c 2 * 2a ~ 2 (5.32) 

In calculating the spread target interference, we made the assumption that the 
refelectivlty of the distributed target in the integral domain is constant. 

5.4.2 AN APPLICATION 

a) In the waveform discussed in the previous section, each pulse in a pulse- 
pair is transmitted with the same carrier frequency, same bandwidth, and same 
envelope but opposite slopes of linear FM. Because the two pulses are transmitted 
with the same carrier frequency, the coherency of the returned signals from the two 
pulses in a pulse-pair is guaranteed if these signals are returned from the same 
resolution volume, The covariance method can be used to estimate the mean 
frequency from each pulse-pair. As discussed in chapter 4. the covariance method 
does not depend upon the inter-spacing between two adjacent pulse pairs. 


T 2 + 5 

C- s 

T 2 -5 

Va-ja 
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Since the maximum unambiguous frequency only depends upon the length 
of T 2 and the maximum unambiguous range only depends upon the length of Tj + 
T 2 , the lengths of Ti and T 2 can be adjusted to achieve a desired unambiguous 
frequency and range. Although the plot of the ambiguity function shows that this 
scheme cannot remove the peaks of the frequency ambiguities, we can still use this 
scheme to get the correct mean frequency of the spectrum. As the two pulses In one 
pulse pair are In the same bandwidth and have the same carrier frequency, the 
echoes from the two pulses can be coherent if the transmitted pulses are coherent. 
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Figure 5.11. Simulation of waveform modulation method for 
estimating first moment with extended frequency range, with one 
run, a=0.3, SNR=10 dB. 


In summary, the waveform modulation allows us to achieve high signal to 
interference ratios and to obtain coherence between the pulses. The radar echoes 
from the two pulses in a pulse-pair are separated through two matched filters. The 
frequency ambiguity is reduced by the pulse pair method with two PRT s (Ti and T 2 ) . 
A computer simulation of this algorithm is presented in Figure 5.11 (with 
measurements as a function of mean frequency). This figure shows that the 
waveform modulation method produces smaller errors in estimation than the 
STAGGER-B method. 
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STAGGER-B method. 


b) Furthermore, we can combine a) with the Inter-pulse coding method. For 
example, we can code the pulse train with a random binary phase coding: each pulse 
pair is multiplied by a random binary phase either 0 or it. One such an example is 
given in Figure 5. 12. By comparing Figure 5. 12 and Figure 5.9, we can see that the 
peaks of the range ambiguities are further reduced by the random phase coding. 
More important is the fact that the spectrum of the interference is flattened. This 
can reduce bias in our estimate of the mean frequency. However, one needs to notice 
that the effective signal-to-interference ratio is not being reduced by the random 
phase coding, only the spectrum of the interference is flattened. 



5.5 RESULTS AND CONCLUSIONS 
5.5.1 THE SIMULATOR 

Two simulators were used in this chapter. One simulator written in C 
language was used to study the ambiguity functions for the dfferent algorithms 
discussed This simulator allows the user to specify parameters such as the pulse 
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length, pulse repetition time (PKR, type of inter-pulse coding, and number of pulses 
in the pulse train, etc. The resulting ambiguity functions are plotted with 3-D 
hidden-line drawings on the computer screen. The second simulator is similar to 
the one used in Chapter 4. This simulator Is used for Monte Carlo simulations of the 
algorithms discussed in this chapter, with specified parameters such as SNR, 
number of runs, spectral width of the signals, etc. 

5.5.2 MONTE CARLO SIMULATION 

Monte Carlo simulations have been used to compare the performance of the 
three algorithms discussed in this chapter: STAGGERED PRF-A, STAGGERED PRF- 
B and waveform modulation methods. As in Chapter 4, the results are presented 
with RMS errors as functions of SNR, spectral widths of the signals, and mean 
frequencies of the signals. 

The data used in the simulations were generated with Gaussian shaped 
spectral density functions with different SNRs and spectral widths. In the 
simulations of the staggered PRF-A method, two different PRF’s were used with a 
ratio of Fi:F 2=7:8. In the simulations of staggered PRF-B method, two PRT’s were 
chosen with the ratio of T2 :Ti= 7:8. In the simulations of the waveform-modulation 
method, the ratio of inter-spacing of the two pulses in a pulse pair and the spacing 
between pulse-pairs is chosen to be 1:7. All the three estimators should have 
unambiguous frequency ranges as large as seven times the Nyquist interval. 

All the computer simulations are based on the normalized Nyquist 
frequency interval defined as [-1.1]. The mean frequencies of the input data change 
from -5 to 5 with specified steps. This frequency range is 5 times the Nyquist 
frequency. Figures 5.14 to 5.16 illustrate the Monte Carlo simulations of the RMS 
errors as functions of the SNR and mean frequency. Figures 5.17 to 5.19 show the 
Monte Carlo simulations of the RMS errors as functions of spectral widths of the 
input data and SNRs. From these results, we can draw the conclusion that the 
STAGGERED PRF-A method produces the smallest RMS errors, followed by the 
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waveform-modulation method. The STAGGERED PRF-B method produces RMS 
errors which are approximately an order of magnitude larger than those of the 
STAGGERED PRF-A method. 

However, the STAGGERED PRF-A method is not without its disadvantages. 
When the spectral width of the signal is large compared with the Nyquist frequency 
(a/ % > 0.3), the STAGGERED PRF-A method produces a significant number of false 
correlations as shown in Figure 5. 16a. The RMS errors produced by the waveform- 
modulation method do not change very much as the spectral width of the data 
changes. 


- 141 - 



PDF of False Estimate 


5 dB SNR 
10 dB SNR 


LO o 


^ LO CD 

J* CD O 


Mean 

Frequency (Hz) 
a) 


OdBSNR 

► 

5 dB SNR 

i 

10 dB SNR 


g g s 

uj ^ o 


, * > a * t ■* i , j t 4 


-6.0 -4.0 -2.0 0.0 2.0 4.0 6.0 

-10500-7000 -3500 0 3500 7000 10500 

Mean Frequency (Hz) 

b) 

Figure 5.13. Monte Carlo simulation of STAGGERED PRF-A method 
with 100 runs, a=0.3, number of pulse pairs = 128: a) probability of 
false estimate; b) RMS error as a function of frequency and SNR . 
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b) 

Figure 5.14. Monte Carlo simulation of STAGGERED PRF-B method 
with 100 runs. a=0.3, number of pulse pairs = 128; a) ensemble 
average of estimates; b) RMS error as a function of frequency and 
SNR. 
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Figure 5.15. Monte Carlo simulation of WAVEFORM MODULATION 
method with 100 runs. ct= 0.3, number of pulse pairs = 128; a) 
ensemble average of estimates; b) RMS error as a function of 
frequency and SNR . 
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Figure 5.16. Monte Carlo simulation of STAGGERED PRF-A method 
with RMS error and probability of false estimate as function of 
spectral width of the signal and SNR, 100 runs, number of pulse pairs 
= 128; a) plot of the probability of false estimate; b) plot of RMS 
errors. 
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Figure 5.17. Monte Carlo simulation of STAGGERED PRF-B method 
with 100 runs, number of pulse pairs = 128; a) plot of the ensemble 
average of estimates; b) plot of RMS errors as a function of spectral 
width of input signal. 
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Figure 5.18. Monte Carlo simulation of WAVEFORM MODULATION 
method with 100 runs, number of pulse pairs = 128; a) plot of the 
ensemble average of estimates; b) plot of RMS errors as a function of 
spectral width of input signal. 
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5.5.3 CONCLUSION 


Among the algorithms discussed in this chapter, the staggered PRF-A 
produces the smallest RMS errors of these three algorithms when the spectral width 
is less than 0.3. The staggered PRF-B method produces the largest RMS errors 
among these algorithms. For all these algorithms, the RMS errors increased as the 
SNRs decreased. For the staggered PRF-A and PRF-B methods, the RMS errors also 
increased as the spectral widths Increased. In the computer simulation, for the 
waveform modulation method, the RMS errors did not change very much as the the 
spectral widths of the radar signals increased. 

In the RAWS system, the spectral widths of the radar echoes are about half of 
the Nyquist interval. In such a case, the PRF-A method has very large probability of 
false measurements, about 20% at 10 dB SNR and 30% at 5 dB SNR The RMS errors 
produced by the waveform modulation are about 7% of the Nyquist interval at 10 dB 
SNR and about 12% at 5 dB SNR The staggered PRF-B method generates too large 
RMS errors to be used In this system. Assume that the PRF is 3500 Hz. Then the 
Nyquist interval is 1750 Hz. Therefore, the waveform modulation method will 
generate about 122 Hz RMS errors in the estimates of mean frequencies at 10 dB 
SNR and about 2 10 Hz at 5 dB SNR 
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APPENDIX 5.A AMBIGUITY FUNCTION FOR WAVEFORM MODULATION 


To derive the ambiguity function for the waveform modulation method, we 
first need to review some basic formulas often used in deriving ambiguity functions 
and the ambiguity functions for some simple waveforms and pulse trains. 

BASIC FORMULAS 


+oo 



-OO 


F( e' st2 ) 


■Vf*- 


0)2 /4 s 


(5.34) 


+oo 


-(at + 2(3t) 



(5.35) 


+oo 

J e ±j2nxa(x-b) dx = ^ b)) (5.36) 

.00 
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For a rectangular pulse, 



The ambiguity function for a rectangular pulse is given in the following form: 
+oo 

X u (*.<l>) = J U W u ( t+x > e dt 

oo 


+ 5/2 

= fu(t) u (t+x) e dt 
. 5/2 



For an FM rectangular pulse. 


uttl^rectQej ™ 1 
The ambiguity function is 
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X u = 


rec Y t\ /(5- It 1 )sin(7r(ctx+<»(S— I x l)) \ 
'28/ ^ rc8(5- 1 1 1) (aT+<j>) J 


AMBIGUITY FUNCTION FOR SI NGLE FM GAUSSIAN PITT .SR 
For a Gaussian pulse 

u(t) =(2a) 1/4 e Kat2 + i Kal2 


The ambiguity function Is 
z u (t,« = e>" x< f e W2 ( a ' 2 + 


AMBIGUITY FUNCTION FOR A PUT.SE TRAIN 
Assume 
N-l 

y(t) = ^u(t-nT 0 ) (5.41) 

n=0 

The ambiguity function for a pulse train is defined as 
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J n=0 
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Assume T q > 25 , if | x - mT Q | 
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(5.38) 


(5.39) 


(5.40) 


(5.42) 
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Xy(M>) 
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(5.43) 
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+oo 


X N-1 

| u(t - nT 0 ) u*(t - (n - m) T 0 + x) e'j 2 * 4 * dt 

n=0 
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Xu ( t + m V« Ze^ 

n=0 
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( 5 . 44 ) 


So, the ambiguity function for a pulse train is 


r n-i 
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AMBIGUnY FUNCTION FOR DUAL WAVEFORM MODULATED PULSE TRAIN 

Assume that y(t) is the sum of two different waveforms u(t) and v(t) where v(t) 
lags u(t) by T. The ambiguity function of y(t) can be written as: 

+oo 

Z y (T.4>) = J(u(0 + v(t-T)) (u*(t+x) + V‘(t - T + x)) e 'J 27t4>t dt 
J _00 
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= X U (T,4» + X^t-T,*) + x vu (* + T -4>) e' i2 ^ T + X vv (x,<t» e^ T 


(5.46) 


Similarly, we can derive the ambiguity function of a pulse pair train. Assume that 
the pulse train can be represented as: 


N-l 

y(t) = ^u(t - nT t ) + v(t - nT r T 2 ) (5.47) 

n=0 

where T 2 < T l and twice the pulse duration 25 <T 2 . Thus, the ambiguity function is 
defined as 
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N-l N-l 


+ I I 

n=0 k=0 


+oo 

J v(t-nT 1 )v’(t-kT 1 -T 2 + t) e') 27t<J)t dt 

OO 


= X £u (t,4>) + X Zv (t,4» ei 2K * T 2 + x Zuv (x - T 2 ,<D) 

+ *Ivu (t + T 2'® e'J 27C<)>T 2 (5.48) 

In special cases, the pulses in a pair are linear FM modulated with an 
opposite frequency modulation as shown below. 

u(t) = (2a) 1/4 e- rat2 + S rtai2 


v(t) =(2a) 1/4 e-* at2 -i* at2 


(5.49) 


Since v(t) = u (t) and from the properties of the ambiguity function, we know that [83) 


X u (x,4>)= X v (x, -40 


(5.50) 


so from (5.40), we can derive Xy(x,<t>). 


X v (x,4>) 


,<)>)= e i rt ^e- n/2 ( ax2+1 ^l 1L ) 


(5.51) 


The cross ambiguity x uv (x,<)>) can be derived from the following integral: 
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+oo 


X^Cc,*) = J u(t) v*(t + x)ei 2n * x dt 

oo 


+oo 

= J(2a)| e'^ 3 " ) a)t2 e' rt(a ' i a)(t + t)2 e ' 2n<|)t dt 

OO 


+oo 

= (2a)2 J e|-2n(a-ja)(t 2 +tT+y)Je'j 27t<|)t dt 

.OO 


+oo 

= (2a)2 e' 7 ^ 3 ' a)t2 J e 2nU " j a)(t2+(x + ^ >‘> dt 

.OO 


If we let a = 2rc(a-ja) and p = ii(a-ja)(x + ), by equation (5.35), the above integral Is 

a-jcx 

equal to 



a-ja 


e *Pl 


,-v, . w 2 2frt(a-ja) -<l> , 

(f<a- a)(-x - z ,) 

2 a 2 + a 2 (a-ja) 2 


"rt. f 2 2xa<t> 2 (a 2 -a 2 )<t> 2 
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(5.52) 


As u(t) 


j2(t>ax j2aa<}> 2 
+ a 2 + a 2 ' (a 2 +a 2 > 2 J. 


■4 


a - ja 


n , 2 

e'2 (aT + 


2 > e'jf (-at 2 ' 2 ^ + ~ 2 * ~ 


a^ + a 


a + a 


* * 

V (t). it can be shown that = x^M). It follows that 


Zvu^.^) 




k / 2 a0 

e '2 (ax + 2 2 
ja a + a 



+ 2({>t + - y^~ ) 
a z + a z 


(5.53) 
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The total interference at the output of a matched filter can be approximately 
calculated as: 


+§+oo 





Aa-iaVt 2 -\ / 2(a 2 + a 2 ) 

* ] y-77fr- e2 a + J a 


^ Je'^-J^dx 
.5 



( 5 . 54 ) 
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Chapter 6 


Conclusions and Recommendations for 
Future Work 

The major focus of this research was a system study of a Doppler radar used 
for global wind measurements. In addition to studying cloud modeling and the 
system configuration, such as antenna scan pattern, we discussed algorithms for 
estimating the first moment of power spectra of radar signals, and algorithms for 
reducing frequency and range ambiguities of Doppler radar systems. The results can 
be summarized as follows: 

(1) In Chapter 2 we reviewed volume backscattering of radar echoes from 
clouds. For three different types of cloud, we simulated the SNRs of the radar 
echoes. The results demonstrated that, from the SNR point of view, frequencies of 35 
GHz or higher are needed to obtain high enough SNRs of radar echoes from these 
types of clouds, presuming the radar system has enough power and enough gain of 
the antenna. Although the results were based only on water-cloud models, ice 
clouds tend to have larger reflectivity than water clouds. The SNRs from Ice clouds 
will be larger than those from the water clouds with the same system parameters. 
However, in the computer simulations, the cloud models were based on an 
analytical drop-size distribution formula, and these drop-size distributions may be 
quite different from the drop-size distributions found In practice. Therefore, more 
realistic cloud models may need to be developed in future study. 

(2) In Chapter 3, we discussed the system configuration of the radar wind 
sounder. Three different antenna scan methods were discussed, and the result 
revealed that a combination of uniform and discrete antenna scanning would result 
In acceptable pointing errors and cost. We also discussed a tracking method to 
estimate the Doppler frequency shift caused by satellite motion. The tracking Is 
through a combination of satellite Inertial navigation system and a second-order 
phase-lock-loop. The result of a computer simulation with a step Input function 
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showed that for given antenna pointing errors, and parameters of the phase-lock- 
loop, this method can achieve very small RMS errors (several Hertz) in a steady 
state. However, the simulation was crude since the clutter echoes used as the input 
signal of the tracking system may be much different from the input signal used in 
the simulation. More thorough analysis of the tracking system and more accurate 
models of the tracking system and input signal need to be developed in future study 
of RAWS. In addition, different tracking methods, such as a Kalman filter, may also 
need to be considered In the future study. 

(3) We also conducted an error analysis in Chapter 3, and the result showed 
that error caused by satellite pointing angle is an important factor to the system 
performance. To achieve the required accuracy (lm s' 1 ), a large number of 
independent measurements needs to be averaged. The error analysis was based upon 
the error In measurement of the mean Doppler frequency and estimate of the 
Doppler shift caused by satellite motion. However, these estimate errors and 
measurement errors are usually functions of the SNR of the received signal. An 
error analysis concerning SNRs of the system may need to be performed In the 
future study. 

We also pointed out in Chapter 3 that clutter rejection is an important issue 
in the design of RAWS. However, in this dissertation we did not discuss this 
problem. This topic should be left as a topic of future study. 

(4) In Chapter 4, we discussed several algorithms for estimating the mean 

Doppler frequency: the FFT estimator, the covariance estimator, and the 

estimators based on autoregressive models. The covariance estimator produced 
slightly smaller RMS errors than the FFT estimator and the autoregressive-model- 
based estimators in the computer simulations. However, these results were derived 
under the condition that the power spectrum of the radar echo is symmetric and has 
only one peak. In practice, interference caused by leakage from the transmitter, or 
clutter, may cause the returned signals to have non-symmetric spectrum or 
spectrum with more than one peak. 
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The performance of the FFT estimator and the AR-based estimators are 
affected by the noise and alias of the frequency of the radar echoes. A noise- 
suppression method was discussed in this chapter. It can effectively remove the bias 
caused by white noise. We also discussed a method for removing the bias caused by 
frequency alias; this method can also reduce the bias caused by noise. 

The FFT, the covariance estimator, and the AR estimator are comparable in 
terms of estimate errors. Although, the RMS errors in the estimates depend on both 
the spectral widths and SNRs of the input signals. The RMS error is mainly a 
function of the spectral width of the signal when the SNR is above 5 dB. As 
discussed, the spectral widths of the signals in RAWS are about half of the Nyquist 
interval. In such case, the RMS errors produced by these estimators discussed in 
Chapter 4 were as large as 10% of the Nyquist frequency. That is, for a 3500-Hz PRF, 
the RMS errors are around 175 Hz. This error is higher than the error limit (23 Hz) 
for the mean frequency estimates needed to produce the lm s' 1 accuracy in wind 
estimates. About 64 independent measurements may need to be averaged to achieve 
the required accuracy. 

(6) In Chapter 5. we reviewed several algorithms for reducing the radar 
ambiguities, such as a phase-coding method and staggered prf methods. However, 
these methods do not perform well when the spectral widths of the radar echoes are 
large ( > 50% of the Nyquist interval). Therefore, a new method based on using 
different waveform modulations on successive transmitted pulses was developed to 
reduce the radar ambiguities. Monte-Carlo simulations were used to compare the 
performance of this method with two staggered-prf methods. The results showed 
that when the spectral width of the signal increases, the RMS errors increase 
rapidly with the staggered-prf methods. However, the RMS error for the waveform 
modulation method does not change very much as the the spectral width of the 
radar signal increases. 

In the RAWS system, the spectral widths of the radar echoes are about half of 
the Nyquist interval. In such a case, the staggered PRF-A method has very large 
probability of false measurements; about 20% at 10 dB SNR and 30% at 5 dB SNR. 


- 161 - 


The RMS errors produced by the waveform modulation are about 7% of the Nyquist 
Interval at 10 dB SNR, and about 12% at 5 dB SNR The staggered PRF-B method 
generates very large RMS errors, over 50% of the Nyquist interval, when the 
spectral width is equal to or greater than half of the Nyquist interval. In all these 
computer simulations, the extended maximum unambiguous frequency is 7 times as 
large as the Nyquist Interval. For a 3500-Hz prf. the maximum unambiguous 
frequency is 12250 Hz. In summary, the waveform modulation method is the most 
promising algorithm among the algorithms discussed in this chapter to be used in 
the RAWS for estimating the Doppler mean frequency. 
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